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FOREWORD 


The objective of this research was the development of an analytical 
technique to evaluate the probability that particles would relocate from 
nonsterile to sterile areas on a spacecraft. This recontamination process 
is important for all multiple missions with separate microbiological burden 
allocations for various major spacecraft systems, and critical for life 
detection experiments that risk contamination from nonsterile components. 

The approach has been to study the effects of typical mission environ- 
ments on the redistribution of particles on spacecraft surfaces both 
analytically and experimentally. This study consisted of three logical 
components, which have been reflected in the effort: (1) particle adhesion, 

(Z) dynamic release mechanisms, and (3) particle transport. The effort in 
particle adhesion has been principally a particle release experimental 
program, together with analytical work and attempts to correlate other data 
found in the literature and elsewhere. Under dynamic release mechanisms, 
meteoroid impact and pyrotechnic firing have been modeled. The particle 
transport activity was an analytical effort which included the development of 
computer codes for spacecraft geometry and orientation, forces acting on 
released particles, and trajectory. 

An analysis of the particle adhesion experimental data was consistent 
with a dominant contribution by Van der Waals (molecular) forces in vacuum. 
These forces have the same linear dependence on particle size as the major 
adhesive force in air with a relative humidity in excess of 65% (capillary 
or water surface tension force). Surprisingly, our numerical results for 
vacuum adhesion are commensurate with published data for experiments 
conducted in air. 

Meteoroid impact has been modeled by the elastic response of a large 
plate to a Gaussian spatial distribution of pressure. The analytical solution 
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for the Green's function is known. The relationship between meteoroid and 
target physical parameters and the time dependence and magnitude of the 
pressure function was developed by an adaptation of an existing analysis. 
Predictions of particulate release due to meteoroid impact were calculated 
and compared with the results of experimental simulations. 

The particle transport analysis was to solve the equation of motion of 
the released particle in the vicinity of the spacecraft. Of the important 
forces, the electrostatic force on the particle due to the electric field of the 
spacecraft is most difficult to predict. Results for the charging rate and 
equilibrium potential of particulates in interplanetary space and an approxi- 
mate electric field were obtained. 

Finally all of these components were assembled into an operational, 
integrated computer code. For a demonstration calculation with this com- 
puter code, a geometrical model of a dual purpose spacecraft and the 
spaceflight phase between Earth orbit and Mars encounter were chosen. The 
results indicate that particulate recontamination is a likely process for 
this mission model. Other predictions, such as the distribution of escapes 
and particulates which relocate on the relatively contaminated regions of the 
spacecraft were also obtained. 
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ABSTRACT 


A computer simulation program to model the redistribution of 
particulate contaminants on a spacecraft after launch has been developed. 
The component models for particulate adhesion, meteoroid impact, and 
electrostatic forces are described and intermediate results are presented. 
The results of a sample calculation have shown that the recontamination 
process is important. 
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I. INTRODUCTION 


Planetary quarantine is an international cooperative program concerned 
with the prevention of the contamination of another planet by terrestrial organ- 
isms. Such a contamination by an automated spacecraft could yield false 
positive results from its own life detection instruments and, if the terrestrial 
microbes grew and spread on the planet, would confuse all subsequent studies 
(Ref. 1). 

In the implementation of a planetary quarantine program for a partic- 
ular spaceflight mission, a distinction is made between spacecraft which are 
intended to probe or to land on another planet and those which are designed 
to flyby or orbit outside its atmosphere (Ref. 2). Thus for a multiple purpose 
mission to another planet, the allowed (and actual) microbial contamination 
of the non-incursive spacecraft components would be significantly larger than 
that of the others. The recontamination process is the redistribution of par- 
ticulates, presumably as sociated with microbes, from a region on a space- 
craft with a relatively large burden to a region that has been decontaminated. 
The objective of the study discussed in this paper is the development of the 
techniques to evaluate the probability of this process. 

The approach has been to study the effects of typical mission environ- 
ments on the redistribution of particles on spacecraft surfaces, both analyt- 
ically and experimentally. This study consists of three logical components, 
which have been reflected in the effort: 1) particle adhesion, 2) dynamic 

release mechanisms, and 3) particle transport. The effort in particle 
adhesion has been principally a particle release experiment, together with 
analytical work and attempts to correlate other data found in the literature 
and elsewhere. Under dynamic release mechanisms, meteoroid impact and 
pyrotechnic firing have been modeled. The particle transport activity was an 
analytical effort which included the development of codes for spacecraft 
geometry and orientation, forces acting on released particles, and trajec- 
tory analysis. 


This study does not treat directly the transfer of viable microbes because 
the relationship is unknown and because of the lethality of natural space 
environments for microbes. 
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Finally, all of these components were assembled into an operational, 
integrated computer code. For a demonstration calculation with this com- 
puter code, a geometrical model was based on a hypothetical spacecraft and 
the spaceflight phase between Earth orbit and Mars encounter was chosen. 
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II. MODELS 


A. PARTICLE ADHESION 

The four types of adhesion forces which may play an important role in 
particulate adhesion to a surface are: Van der Waals, contact potential, 
Coulomb, and capillary. The Van der Waals or molecular force is propor- 
tional to the particle dimension, depends critically on surface roughness, 
and is fairly unaffected by material conditions. Although the other forces 
may be much larger under certain circumstances, the molecular force is the 
major source of adhesion for inert contaminants seeded in air onto a surface 
and then evacuated. Contact potential differences due to the surface effects 
of dissimilar materials cause an electrostatic binding proportional to the 
particle dimension to the two-thirds power. Coulomb forces arise from 
actual charges on the particles due to external ambient electric fields. Since 
these forces are inversely proportional to the particle dimension, they are 
quite important for small particles. The electric charges involved, however, 
in both contact potential and Coulomb forces will neutralize in the presence 
of ambient water. More importantly, for particles seeded in air (or even in 
dry nitrogen) onto a surface which is then evacuated, the Paschen limit 
places an upper limit on the residual charge. During evacuation, a region 
of corona breakdown, where the field in the small gap between particles and 
the surface causes the air to ionize, is passed. Finally, the capillary or 
water surface tension force, which is also proportional to particle size, is 
potentially largest. It is, however, negligible in a reasonable vacuum. At 
ambient air pressures, the capillary force exhibits a definite hysteresis 
effect with respect to relative humidity. That is, the force depends on the 
history as well as the relative humidity at a given time. 

Some conclusions relevant to the recontamination task may be drawn 
from the preceding discussion. The Van der Waals force, proportional to 
particle size, is the best model for the vacuum problem in space. During 
the ascent, the adhesion force changes in a very complicated fashion, but it 
approaches the Van der Waals as a lower limit. Conversely, a particle 
adhesion experiment in vacuum should reliably measure the Van der Waals 
force. One may compare vacuum results with results obtained under 
ambient pressures only with caution. 
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Previous models used in this study have related the applied force (given 
by the product of the applied acceleration and the particle mass) to the 
removal fraction. In each case a characteristic acceleration or force corre- 
sponding to a fixed removal fraction was defined. In Ref. 3, the character- 
istic acceleration a Q for a removal fraction of 0.63 was given as: 


2 X 1Q" 4 (0,4 + 0.006 RH) 
TTd^pG 


( 1 ) 


2 -3 

where a Q is in units of kilo-gee (kG), p is the particle density (g pm ), d 

is the particle diameter (pm), G is the acceleration due to gravity 

(980 cm s“ 2 ), and RH is the relative humidity (percent) during release. 

The numerical constants in Eq. 1 were determined empirically from data 

in the literature. 

Subsequently, an experimental program was conducted to obtain data on 
particulate adhesion under vacuum conditions. The' data to be modeled 
(Figs. 1-5) were obtained by an impulse method for glass beads on stainless 
steel. The details of the apparatus and the test procedures have been previ- 
ously reported (Refs. 3 and 4). It should be noted that in these tests the 
glass beads were seeded in air, and then the removal occurred in vacuum. 

As discussed above, for this experiment and for the recontamination 
analysis, one expects a characteristic force, F corresponding to a removal 
fraction of 0.5 given by: 


F = kd 
o 


or the characteristic acceleration a 

o 


given by: 


a 

o 


6k 

,2 

Trd p 


(2) 


(3) 


2 2 

^One gee equals 9.8 m/s , the acceleration due to gravity at the surface of 
the earth. 
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where d is the diameter and p is the density of the particle. Note that 
Eqs. 1 and 3 are of the same form for zero relative humidity. 

The next component of the model must relate the applied force or accel- 
eration and the characteristic force or acceleration to predict the removal 
fraction. Since a distribution of adhesion forces is expected even for a 
collection of identical particles because of the variation in microscopic sur- 
face conditions, a probabilistic model is indicated. In the present model, the 
constant k of Eq. 2 is assumed to be log-normally distributed with a mean 
value of log k equal to m and a standard deviation cr^. Then the removal 
fraction for a specific d is given by the probability that the applied force F* 
exceeds the characteristic force F . This probability may be written: 


P(F* > F q ) 



(4) 


where 


= F*/d = 


,2 < 

ird pa 


(5) 


and a is the applied acceleration. 

By inverting Eq. 4 for a given data point and an assumed value of ov 

one can calculate an estimate for m. Such a collection of estimates for m 

may be averaged to provide a best-estimate m Q for the selected u* . There 

is also a standard deviation of this collection, a* , which should be mini- 

mized. The program for search over cr , treating the data, and minimizing 

o' has been incorporated into a computer code. The best fit to the data was 

obtained for c r = 0.4, m = -0.882, and cr = 0.52. These values are 
r 9 o 9 m 

consistent with MKS metric units, e.g., k"‘' and 10 m ° have units of newton 
meter ^ . 

Up to this point the theoretical development is due to B. Nelson (Ref. 5). 

Nelson, however, uses Eq. 4 with cr substituted for cr and m = m for 

the most probable removal fraction. The correct procedure appears to be 

the use of cr and m for the most probable calculation. Further, the 
r o * 
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meaning of <r is such that if one calculates the removal fraction with 

m = m - tr substituted into Eq. 4 and with m = m + cr , there is the 
o m o m’ 

usual 68% probability that a test result will fall between the two answers. If 
one is concerned only with an upper limit on the removal fraction, there is 
an 84% (50 + 34) probability that the prediction with m = m Q - <r will not 
be exceeded. These comments may be extended in the usual way. 

The results of this analysis are shown in Figs. 1-5. Each figure cor- 
responds to a particular size of glass bead. The data and the predictions of 
Eq. 4 are given as removal fractions in percent vs. surface acceleration in 
kG. The most probable and the "conservative" (m = m Q - <r ) predictions 
are shown and, in the case of 110-pm beads, the "non-conservative" 

(m = m Q + cr m ) prediction is also given. 

It is interesting to compare our new results for the adhesion force with 
the previous model and with results published in the literature. The present 
prediction of the mean adhesion force (N) is given by Eq. 2 with the diameter 
d expressed in meters and a value of k = 0.13 N/m. Corn (Ref. 6) used the 
same formula and obtained the following values (converted to MKS units): 

0.17 N/m for the adhesion of Pyrex particles to optical glass, 0.12 N/m for 
the adhesion of quartz particles to optical glass, and 0.075 N/m for the 
adhesion of quartz particles to Pyrex glass. Since these data were obtained 
in air at 95% relative humidity, the agreement with the present results must 
be considered somewhat fortuitous. 

A comparison with the characteristic acceleration of the previous 
model, Eq. 1, requires some algebra. One may express the prediction of 
that model in terms of the force required to cause a 63% removal in air 
(RH = 0%): 


F(0.63) = 0.13d (6) 

where F is expressed in newtons and the particle diameter d is expressed 
in meters. The prediction of the present model has a coefficient of 0.18 for 
this removal fraction (compare to 0.13). This fair agreement may also be 
coincidental. The previous model was based on data for many materials in 
air with a relative humidity of 50 - 60*%. 
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From a comparison of the data in the literature and a consideration of 
the Van der Waals and capillary adhesion forces only, one may predict the 
trend of total adhesion between ambient pressure with 100% relative humidity 
and vacuum. As the relative humidity is decreased, the capillary force 
decreases in a linear fashion (Ref. 6). As long as water is present, the Van 
der Waals force is reduced to a negligible factor by the molecular interaction 
with the interface water (Ref. 7). For small relative humidity (<50%), there 
is a minimum adhesion. Finally, as the test system is evacuated, the capil- 
lary force becomes negligible and the Van der Waals force attains its maxi- 
mum value. On the basis of the data, the adhesion force in vacuum is com- 
mensurate with that in air at 100% relative humidity and certainly exceeds 
values for low (but non-zero) relative humidity. 

This particle adhesion model has been developed into a computer code 
called RELEAS. A version modified for its integration into the complete 
recontamination code appears in Appendix C. An independent version which 
accepts surface response data in punched card form is also available. 

B. METEOROID IMPACTED SURFACE RESPONSE 

The meteoroid impact model is intended to provide surface response 
characteristics as a function of meteoroid and surface physical parameters. 

The connection between the meteoroid impact parameters and plate 
response is the loading function. The loading function model predicts a 
pressure P, in general a function of time t, and r, position on the surface 
relative to the center of impact, for a given meteoroid event. For simplic- 
ity, the positional dependence has been assumed to be Gaussian and 
factorable: 


2 . 2 

P(r, t) = P Q e' r /S f(t) (7) 

This choice was motivated by the existence of a closed-form solution 
for the velocity response of a thin plate subjected to a loading of this form 
with the impulse time function (Ref. 4). This response function is related 
to the response due to a general loading time history f(t) by convolution 
integrals. 
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In the preliminary considerations, it was found that for reasonable 
forms of the time dependence f(t), only a characteristic time T q of the 
impact was important for the peak acceleration, our primary interest. The 
detailed behavior of f(t) certainly determines the detailed time history of 
the surface acceleration. 

The function f(t) has been chosen: 


f(t) = 



( 8 ) 


This form has been derived along the lines of an analysis due to Ludloff 
(Ref. 8). This treatment takes into account thfe transport of molten target 
material in the crater region. The result depends on Ludloff's form for the 
crater radius; the crater is assumed to be hemispherical. The expression 
for the crater radius R in terms of the target material strength S, the pro- 
jectile diameter d, density Pp and velocity v q is: " 


R 1 
d “ 2 



( 9 ) 


This formula agrees reasonably well with the data in the literature. 

In this analysis, the time history of Eq. 8 is valid only for the boundary 
of the crater being formed. This restriction appears to pose no problem in 
calculating the post-impact elastic wave far from the crater. The solution 
of the equation of motion leads to Eq. 8. The characteristic time T q may be 
immediately identified as: 



4p, 


( 10 ) 


One may also identify P q a 2S. In the formula for T p T is the target 
density and all other symbols have been previously defined. Ludloff, in 
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keeping with his treatment of the molten target material, takes as the target 
material strength S: 


S = Prj, 


(ID 


In this formula X.^, is the latent heat of fusion. For aluminum, for example, 
S = 8.7 X 10 8 N/m 2 . 

Finally, we have chosen for the parameter s in the loading function 
(Eq. 7), a value of (1/3)R, where R is the crater radius given in Eq. 9. At 
this point the loading function is completely expressed in terms of known 
parameters of a meteoroid impact. 

As stated previously, classical plate theory provides the velocity 
response WDOT of a thin plate to an impulsive loading with a Gaussian shape 
factor. A computer code written by J. Yang yielded the surface velocity w 
and surface acceleration w by a convolution with f(t): 


w(r, t) = / dx WDOT(r,T) f(t - t) 

o 


w(r, t) = / dx WDOT(r,x) f(t - x) 

o 


( 12 ) 


(13) 


In the process of inserting f(t) into this thin-plate code THINPL and 
testing it, a close inspection of the code was made to find a way of shorten- 
ing it. By an analysis of the convolution integral and the shapes of WDOT 
and f(t), we have succeeded in predicting the time t for a given lateral dis- 
tance r at which peak acceleration occurs. Since the calculation of the sur- 
face velocity and acceleration for a given position requires two integrations 
for each value of time, a great savings is realized by limiting the values of 
time to near the time of peak acceleration. In addition, the convolution inte- 
gration can be optimized for these values. As a result, the present code 
yields peak acceleration and the simultaneous velocity as a function of posi- 
tion by a limited search over t. An entire curve of peak acceleration versus 
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position is computed in about the time the previous code required for the 
time history of the acceleration at one position. 

The analysis of the impact of a small meteoroid and/or the analysis for 
positions close to, but not in, the impact area entails a consideration of the 
semi-infinite half-space problem. This analysis, which treats the plate as 
though it were infinitely thick, has been published in the literature (Ref. 9). 
The computer code based on this analysis yields peak accelerations inversely 
proportional to the position r beyond the impact locus. Because our prob- 
lem is limited to finite plate thicknesses, the absolute values of the peak 
accelerations predicted are suspect. 

In the past we have taken the transition between the "thick" and "thin" 
regions to be the value of r equal to three plate thicknesses and have matched 
the peak acceleration at that point. In fact, however, this transition region 
beyond which shear waves dominate the compressional waves depends on the 
other plate properties and the impacting projectile properties. Fortunately, 
the prediction of the classical thin- plate program levels off and then actually 
decreases with decreasing values of r. One may interpret this odd behavior 
of the peak acceleration as due to destructive interference between the var- 
ious contributions to it from the different parts of the finite region of dis- 
turbance. Since the thin-plate analysis treats shear waves only, the leveling 
off of its prediction identifies the transition region uniquely. 

The meteoroid impact program developed for this task, YANG1/ 
THINPL, uses this identification by searching for the maximum peak accel- 
eration predicted by the classical thin-plate program for decreasing values 
of r. For all smaller values of r than the specified value found, the peak 
acceleration is taken as inversely proportional to r. The infinitely thick- 
plate analysis code, a cumbersome long-running code, is not used at all for 
production runs. For very large values of r, the prediction of the program 
drops off finally below a cut-off in acceleration, and calculations cease. In 
addition to the version of the meteoroid impact code YANG1 /THINPL used in 
the complete recontamination computer program (Appendix C), an independent 
version which produces the surface response summary in printed, plotted, 
and punched card form is available. This form of YANG1 /THINPL may be 
used in conjunction with the independent particle adhesion code RE LEAS to 
provide predictions of particle removal due to meteoroid impact. Such 
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predictions, shown in Figs. 6-9, compare favorably with the experimental 
results of a group at Langley Research C enter (LRC) (Refs. 10 and 11). A 
more complete discussion is available in R,ef. 3. 

C. ELECTROSTATICS 

The particle transport analysis addresses the problem of following the 
motion of the released particles. Formally stated, given the particle param- 
eters, the initial conditions, and the forces acting on the particle, the analy- 
sis is to solve the equation of motion. Of the important forces, the elec- 
trostatic force on the particle due to the electric field of the spacecraft is 
most difficult to predict. A major effort within this task has been in the 
electrostatic area. The basic approach has been to analyze particle charging 
so that the particle charge as a function of time is known and to analyze the 
potential solution for the spacecraft- solar plasma system (Ref. 12) to develop 
an approximate electric field. 

The particle -charging model deals with spherical particles for simplic- 
ity and consists of three special cases: illuminated particle, shaded (by the 
spacecraft) particle in the wake, and shaded particle outside the wake. The 
wake region of the spacecraft is a complicated analysis that is discussed 
under the electric field heading. Given the geometry of the particle and the 
potential -dependent currents flowing into it, one may express the charging 
rate in an implicit form. 

For the illuminated particle, the charging rate d4>/dt is given by: 


W = 2 * 9 ^2 10 Ud [ A p + (1 " a) A e + A v] 


(14) 


where 


A p = n p v o /4 for |<!>|« Mv 0 / 2e 


exp (e«t>/kT e ) 


- n v /4 X 
e e 


1 + 2e 


+/kT e ] 


1/2 


for c|> < 0 


for 4> > 0 
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photoelectron flux, dependent on <j> 
secondary emission yield, dependent on <j> 


A = 

v 


a = 


-3 

In this equation, the proton flux depends on the number density n^m ) 

and the directed velocity v (ms *), and the electron flux A depends on the 

-3 ° _ - 1 e 

number density n (m ), the mean thermal speed v (ms ) and the tempera- 
e e 

ture kT e (eV). The photoelectron flux is calculated from photoelectric yield 

data (Ref. 13) for the solar spectrum, material of the particle, and surface 

-2 - 1 

potential 4> (V). With A , A g , and A v in units of m s , the particle 
diameter d in m and the heliocentric distance R in AU, one obtains d<J>/dt 
in Vs -1 . 


A model for the secondary emission coefficient has been devised based 
on an empirical formula for the secondary yielcl of a given material as func- 
tion of the primary electron energy (see Ref. 14). This yield function was 
formally integrated with an assumed Maxwellian solar -wind electron energy 
distribution with a specified temperature. The resulting coefficient is the 
ratio of the secondary electron flux to the primary (incident solar wind) 
electron flux. Parameter values for the empirical formula were also taken 
from Ref. 14. The dependence of the coefficient on the particle potential is 
of the form: 


a- = (A - Be<|>)/(kT e - e<|>) 


for <j> < 0 


(15) 


For the shaded particle, the photo electron flux is zero. 


The electron 


flux A g and the secondary emission yield are essentially unchanged. The 
crux of the problem is the proton flux A^. The proton accretion may be 
expected to depend only on the thermal motion of the protons, since there is 
no line of sight in the shade. It follows then that for this case; 


exp ( -e4>/kTp) 


n v 

-P P X 


2e 


vkT ; 


1/2 


for <J> > 0 


for 4> < 0 


( 16 ) 
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where is the proton mean thermal speed and kT^ is the proton temperature. 

This equation, an analog to the electron flux equation, is only approximate. 

The quality of the approximation depends largely on the value taken for n , 

P 

an effective proton number density in the shade. 

The effective proton density n in the shade depends on the complete 

- P 

solution to the potential problem of the spacecraft itself. Poisson's equation 
relates the potential everywhere to the charge distribution everywhere. The 
boundary condition (surface potential of the spacecraft), the potential, and 
the charge distribution form a self-consistent set. For an equivalent spheri- 
cal spacecraft, a conical wake region is assumed where n =0. The base 

P 

of this cone contains the center of the sphere and is a cross-section of the 
sphere. The altitude of the cone lies in the anti-Sun direction. The half- 
angle of the cone is given by (Ref, 15): 


tan 0 



(17) 


Outside this cone n^ is taken to be n^, the undisturbed proton number density. 


The particulate charging rate for the three conditions (sunlight, shade, 
and wake) is calculated in a computer code PARPOT (Appendix C), as well 
as the necessary terms like the photoelectric current and the secondary 
emission coefficient. An independent version of this code, which provides 
values for these quantities, and also the equilibrium potential of a particle 
in the sun is available. Sample results are shown in Figs. 10 and 11. 


The analysis of the spacecraft electric field or potential problem 
requires many of the same considerations as the particle charging model. 
However, the spacecraft surface potential is merely a boundary condition 
rather than the answer. Fortunately the spacecraft/plasma system equili- 
brates rapidly (time proportional to the inverse of a size dimension) so that 
a time -independent equilibrium analysis is adequate. The problem of the 
potential distribution of an object in a plasma, the plasma sheath problem, 
is known mathematically as a Poisson-Boltzmann or, in the case of a non- 
interacting plasma, a Poisson-Vlasov system. 


The three-dimensional potential solution for an object moving relative 
to a plasma is very difficult even in the simplest geometry, a uniform 
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conducting sphere. This situation, which has been solved implicitly by 
others, requires that the directed relative velocities (as opposed to the 
thermal, random velocities) of all components of the plasma be equal. 

These conditions are satisfied when all of the directed relative velocity is 
due to the object's motion. In the case of interplanetary spaceflight in the 
solar -wind plasma, however, the proton component has a directed velocity 
which dominates its thermal velocity, while the electron component by com- 
parison is essentially at rest. Finally, the photoelectric effect produces yet 
another plasma component exactly at rest with respect to the spacecraft. 

Thus the potential problem is a much more difficult one than a mere compli- 
cated geometry to be approximated by known solutions for simple geometries. 

The general approach taken was to artificially sector the region of 
space about the spacecraft into the near -field, wake, and far -field sub- 
regions (Fig. 12). The electric field in the near-field case was approxi- 
mated by a one -dimensional solution for a plate of material nearest the given 
position. The near -field boundary was taken as the characteristic e-folding 
distance (effective Debye length) for that material. The electric field model 
in the wake regions was adapted from an analysis by Al'pert et al, (Ref. 16). 
The wake boundary was taken to be a cone with a half -angle given by Eq. 17 
(see previous discussion on particle charging) and a base radius (XYWAKE) 
calculated on the basis of a circle with an area equivalent to the spacecraft 
shade projection. The far-field electric field model was an equivalent sphere 
with a uniform potential equal to the area -averaged value PHIAVE and with an 
effective e-folding parameter equal to an area-averaged value ALAMAV. 

The solution to the one -dimensional Poisson- Vlasov problem of a plate 
exposed to solar illumination, electrons, and protons has been found in the 
literature (Refs. 17 and 18). In this solution, the positional dependence of 
the electric potential in equilibrium was expressed implicitly in terms of the 
surface potential and the plasma parameters. An approximate explicit form 
for the potential based on this work has been developed which provides an 
approximate formula for the electric field as well (Ref. 19): 


E = 4c 2 {- Cl + 


-.21 


(jv^-r - °2 x ] 


- C 2 X 


(18) 
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for the region 


x < 


(h Q ' 

C 2 


+ c 


l) 


1/2 


(19) 


where and are known constants, which are calculated. 

In order to obtain numerical values for the electric field, the plasma 
parameters and the equilibrium surface potential are required. The perti- 
nent parameters of the solar-wind plasma are available (Ref. 20). For the 
photo -electron plasma, the typical experimental values obtained for ultra - 
clean surfaces are expected to differ drastically from values corresponding 
to realistic spacecraft surfaces. However, a detailed experiment for real- 
istic surfaces of several spacecraft materials has been reported (Ref. 13). 
Finally, the equilibrium surface potential for these materials has been deter- 
mined with the use of the parameters above and the requirement that the net 
current becomes zero at equilibrium. The values of these parameters are 
given in Table 1. 

Following the authors of Refs. 17 and 18, two classes of solutions to 
the problem were found and noted as monotonic and non -monotonic. The 
quantity of interest in the electric field in the vicinity of the plate is, how- 
ever, virtually the same for the two solutions. Some results for aluminum 
and silica are shown in Fig. 13. The most striking feature of these results 
is a falloff of the field with distance at a far faster rate than the field for the 
case with no photo-electric effect. The latter problem has an exact solution, 
a decreasing exponential, with an e -folding distance equal to the Debye length. 
For the solar -wind plasma parameters, the Debye length is about 12 m. In 
contrast, the e-folding distance, as determined by an exponential fit to the 
results given in Fig. 13, is about 0.8 m for aluminum and 2 m for silica in 
our analysis. Thus, the effect of the photoelectrons is to dramatically 
reduce the effective range of the electric field and to render the effective 
range of the electric field material-dependent. 

The one -dimensional problem in the shade was modeled by a surface 
potential equal to minus three times the electron temperature (in eV) with an 
exponential spatial dependence. The Debye length in the absence of photo- 
electrons was employed. 
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The electric field for the wake is approximate and basically depends on 
a potential model of the form: 


4 > 



( 20 ) 


where (r, z) is the position in cylindrical coordinates , a is the cone radius 
at height z, <j> Q is the shaded surface equilibrium potential, and J q is the 
zeroth order Bessel function. Corrections were made to approximately 
satisfy the boundary conditions of a conical geometry. 

The electric field in the far-field region is based on the usual potential 
solution for a uniform sphere in a plasma (at rest) with a known Debye length: 



( 21 ) 


In the present application, the surface potential <J> o is identified as 
PHIAVE, the sphere*s radius a is identified as XYWAKE, the Debye length 
is identified as ALAMAV, and r is the spherical coordinate distance. The 
electric field may be obtained by a formal differentiation of Eq. 21. 

The calculations of the various surface parameters and averaged 
quantities occur in the computer code ESURF (Appendix C). The determina- 
tion of the appropriate model for a given spatial position and the electric 
field evaluation are accomplished in EFIELD (Appendix C). 
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Table 1. Potential parameters 
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Fig. 1. Comparison of particle adhesion data with predictions of 
model, 22-pm glass beads 
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Fig. 2. Comparison of particle adhesion data with predictions of 
model, 48-p.m glass beads 
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Fig. 3. Comparison of particle adhesion data with predictions of 
model, 68-pm glass beads 



Fig. 4. Comparison of particle adhesion data with predictions of 
model, 94-pm glass beads 
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Fig. 5. Comparison of particle adhesion data with predictions of 
model, llO-pm glass beads 



Fig. 6. Glass bead profile predicted post impact for 
simulation of LRC experiment. Case 1 
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Fig. 7. Glass bead profile predicted post impact for 
simulation of LRC experiment. Case 2 
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Fig. 8. Particulate removal simulation of 
LRC Run #6 (NASA TN-7494) 


JPL Technical Memorandum 33-737 


21 




CHARGING RATE -7-, V/s 


00 


50 GLASS BEADS 


• RECONTAM SIMULATION 
(AVERAGED OVER R) 

□ IRC RUN 6 (FIXED R)* 



Fig. 9. Release velocity distribution simulation of 
LRC Run #6 (NASA TN-7494) 
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Fig. 1Z. General indication of the spacecraft 
electric field with component material designations 
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Fig. 13. One -dimensional electric field 
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III. COMPUTER SIMULATION PHILOSOPHY 


A. BLOCK DIAGRAM OVERVIEW 

The integration of the specific models for the spacecraft surface -impact 
response from micrometeoroids, the particle release prediction, and the 
pyrotechnic device firing is displayed conceptually in Fig. 14. These major 
computer -coding subgroups indicate the general relationship among the con- 
stituent modeling efforts. 

The Meteoroid Model block represents the establishment of a ten by 

6 ) 1 3 

three matrix containing ten mass groupings which span 10 to 10" kilo- 

4 4 4 

grams and three velocity groupings of 1.6 x 10 , 2.4 x 10 and 3.2 x 10 m/s. 
These parameter ranges will later be seen to encompass meteoroids of real- 
istic concern as discussed in Ref. 21. The output matrix contains the come- 
tary meteoroid total mission fluence for each mass-velocity (M-Y) matrix 
location and the cometary meteoroid M-V group flux. 

Once the flux and fluence calculations have been made, the data for each 
M-V group are transferred to the Surface Impact Model code. At this point, 
the effect of a single meteoroid impact upon a thin plate of specified thickness 
is determined with the output formulated into three arrays of data, corre- 
sponding to each M-V sample. These data arrays contain the peak accelera- 
tion of the surface and the velocity of the surface at the time when peak accel- 
eration occurs, each as a function of the range from the impact (third array). 

The simulation has chosen sample meteoroids that span the range of 
physical interest in the case of particulate redistribution. Other meteoroid 
models could also be introduced that yield output in a suitable form: accel- 

eration and velocity of the surface vs. range from impact and a representa- 
tive fluence evaluation. This is essentially what the PYRO EVENT Model 
grouping produces for the relevant pyrotechnic devices of known response 
and location (Appendix B). 

The next group of computer coding represents the prediction of particle 
(on the impacted surface) release probabilities that depend upon the surface 
environments described above. The internal probabilities are functions of 
the particulate diameter (10 size groups) and the range from impact. The 
properties of symmetry and reciprocity allowed all impact ranges to be 
tested at any one position. 
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The number of trial ejecta histories required for full simulation has 
been increasing during each step in the discussion. For each of the thirty 
M-V groups, ten ejecta sizes were analyzed along with, typically, six to nine 
ejecta velocities (see Appendix C). This totals approximately to 2500 trial 
histories necessary to fully cover the variable ranges developed to this point 
in the program flow. 

The next consideration concerns the location of impact on the geomet- 
rically specified spacecraft. Reference 21, upon which the cometary mete- 
oroid information has been based, describes an isotropic directional distri- 
bution; thus, all portions of the exposed spacecraft are assumed equally 
likely to have micrometeoroid impact. The particle transport code, there- 
fore, begins with a Monte Carlo selection of the surface position, limited 
to the unsterile surfaces (Fig. 15), to be considered as the initial trajec- 
tory location. The Monte Carlo method has the advantage of producing 
an unbiased set of random impact points. A disadvantage of Monte Carlo is 
that a large number of selections must be made to obtain statistically repre- 
sentative data. For this reason, Monte Carlo was applied only to this one 
parameter, the position. An analytically deterministic approach using the 
discrete group distributions generated above for all other relevant param- 
eters was applied for statistical and economic optimization. 

The program is arranged such that the descending order of selection 
within the coding follows: M-V groups; number of positions to be considered; 

ejecta size; and velocity groups, respectively. The arrangement allows the 
testing of all ejecta size and velocity groups (60-80 total) at every position 
for each M-V set (2500 times number of positions). 

The final computer code group keeps account of statistical informa- 
tion as the model courses its way through the execution. During each 
kinematical grain trajectory, positional track is kept on the ejected particle 
to monitor for spacecraft recontact and recontact location. The categories 
of interest are whether the ejecta recontacts a nonsterile surface, a sterile 
surface, or escapes. 

B. REQUIRED INPUT DATA 

The spacecraft being modeled is required to be of the basic configura- 
tion displayed earlier in Fig. 12. The truncated cone dimensions are 
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arbitrary, in that the inclination angles and base dimensions are adjusted 
to describe a specific mission. Figure 15 identifies the spacecraft model 
components. The spacecraft surface must be specified for density, thick- 
ness, particulate load, and physical dimensions. Solar generated plasma 
parameters, solar energy spectrum, material photo-electron yield, physi- 
cal constants, mission-dependent information, etc. , need characterization. 
Appendix A lists each specific data entry requirement within the coding for 
element DATAIN. 

C. OUTPUT DATA FORMULATION 

Upon completion of the mission data entry and the compilation of an 
executable element {Appendices A and D, respectively), the execution com- 
mences. The following extended list, Table 2a-i, represents samples of a 
typical computer run printout. Initially, the mission-dependent internally- 
calculated constants are listed as in Table 2a; the tangents of the cone section 
defining angles are listed; the areas of the various sections (i, e. , ASP is the 
solar panel area, etc.); the Z axis cone truncation intercepts are listed, with 
Z1 defining the dark side Debye length; the radius and depth of the plasma 
wake (XYWAKE, ZWAKE, respectively); the Z axis coordinates of the cone 
verticies (AL); the electrical characteristics of the spacecraft, where 
ALAMB is sunlit side Debye length for aluminum and silica, EFEL is the 
sunlit side electric field magnitude, EFED is the shaded side electric field 
magnitude, PHIEL and PHIED are the illuminated and shaded surface poten- 
tials for aluminum and silica, DEB is the shaded side Debye length, ALAMAV 
and PHIAVE are area averaged Debye length and surface potential (for far- 
field considerations), respectively; and finally, the AMAT array describes 
the spacecraft material beginning with the element A top plate and ending 
with the element B bottom plate (see Appendix A). 

The next page in the printout shows the data card image for the case to 
be studied (Table 2b). We have listed the image for a typical meteoroid 
mass -velocity group case study and the image for a pyrotechnic event simu- 
lation (see Appendix E). The two parameters designated as pyrotechnic 
localizers limit the random pyrotechnic position selection to a particular 
section of the space vehicle (see Appendix B). The next page (Table 2c) lists 
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the group identifiers, meteoroid mass, meteoroid velocity, and M-V group 
flux and fluence values. Immediately thereafter, the table of range, peak 
acceleration vs. range, and velocity vs. range for the impacted surface is 
printed. (Table 2c data correspond to points A and B on Fig. 14.) In the 
pyrotechnic event phase, the corresponding page output is shown in Table 2d. 
This likewise corresponds to point B on Fig. 14. 

Preliminary summaries in the form of "Box Scores" are then output 
during execution for each position selected (part D on Fig. 14). Tables 2e 
and 2f depict this information for one position of a meteoroid M-V group and 
a pyrotechnic event, respectively. The position coordinates are for the initial 
ejecta location (Z axis origin is one Debye length in the anti-Sun direction). 
The Box Score indicates the number of ejected grains that recontacted a safe 
area, recontacted a "sterile specified" zone, passed near but missed the 
spacecraft, and escaped at this position. The average velocity of escape in 
m/s is also indicated. 

Table 2g shows the altered form of the output when a re contamination 
event occurs (part C on Fig. 14). This output indicates the ejecta size and 
velocity, as well as several constants relating to the particular character- 
istics of this grain, followed by a step-by-step account of the trajectory from 
meteoroid impact triggered grain ejection to the sterile zone recontact. This 
accounting includes the integrated path length, the positional coordinates, the 
particle charge, the visually helpful shade-solar panel-sector indicies, the 
electric field vector, and the incremental time -of -flight from previous posi- 
tion. Upon recontamination, a recap is printed with a Box Score list for this 
position. When all of the positions have been tested for one meteoroid M-V 
group, another Box Score is printed, as shown in Table 2h, that reflects the 
total mission results due to this particular meteoroid group. 

The general software is configured to study many M-V groups in 
sequence, after which the grand ensemble values for all the groups studied 
are listed. Table 2i concludes the printout showing the total mission fluence 
due to these M-V groups, the grand ensemble, the standard deviation for the 
recontamination prediction, an angular distribution listing for the grains that 
escape, and a probability value (assuming Poisson statistics) for the event 
that one or more grains recontaminate the sterile zone. 
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Fig. 14. Major software component block 
diagram 
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Z AXIS 



Fig. 15. Cone spacecraft model cross-sectional 
geometry 
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TABLE 2a. Initial constants for hypothetical mission (sample computer printout) 


T*( fOtlO«IN4 tut 4t«e»*L SPACECRAFT NODCL OATA 


Thc TaNCCnTS fOR tout AN6LCS ARE *7*44 
Alp«2Q*B ATP*I««2 AT S a l 9 «Q 

Urn l§.22 22* |9*3 r 23* |9*79 

UlAKC* 2MI7V 2 «aK£" 4*9* 0 3 


-*4937 r9*l99 7 • 1 94 

A | *34 • 4 A | J *7 ATI 

29* 20*4* l%» 22 v 4| 


IMC AL ARRAY CONTAIN* j7,9* 20*41 

• • AIaNI • • • • |F(L * • 

• 74J | • 72 1 9 # | 2*3| 


24*01 4*990 

• • £reo • • 

*3*29 *3*29 


• • PHlIl • • 

1 3*0 J * 49 


• • PHUO • • 
•40*0 *40,0 


• AlAMAV • • PH I A v £ 
9*44 *24*2 


OCB 

l «*2 


PR | • | 74 


ABC9«U«7 


At At M02 AC AL AL *|02 S|Q2 S|02 

• * el array • • 

At S 1 02 


TABLE 2b. Typical initialization data printout (sample computer printout) 


The input hetcoro i o case data c*ko for hass group 9 and veloc* group 2 


1 22 1 10 


• 1 0000* 1 0 2 


9 * 4294*09 


OUMMY 

VARIABLI 


IMPACT 
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IMPACT 
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PYRO NUMBER 
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PYRO LOCALIZING 
PARAMETERS 


SECTOR 
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TABLE 2c. Meteoroid model identifiers and surface response model data 

(sample computer printout) 


MASS GROUP* 9 YELOC* GROUP* 2 

the meteoroid ma** I* i*oooo-n the meteoroid veloc* i* 2.i*»**09 

KITH A GROUP Flu* OF 9 a 978*-07 AND THE TOTAL NETEoRqjO IMPACTS uF T H j S TYPE ExPECTtO ARE 1092* 


20 
2 I 
22 


Rll i 

2* 0090**09 
9*o§o9u*Q9 
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3«72l§ f *03 
9*92§73-03 
§*|29»/*03 
§'•3901*03 
**§3*1§*Q3 
7# 2922**03 
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••93**9-03 
1 * 00000-02 
2*00000*02 
3*00000*02 
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*•00000*02 
••00000*02 

1 *00000*01 

2 tOOOOO-0 1 
9*00000*01 


PKACC* I 1 
2*17 l**+Q 1 
9 • *2 l | 2*00 
2 • 7 | 1 G I *00 
I • ••* 1 7 *0Q 

I • 99*02+00 
I • 17293*00 
9 *•••99*01 
••§0G*9-0| 

7*979o*-OI 
*•*73*2*01 
*•02977-01 
§••29*7-01 
S*0l9«t-0I 
9.97 1 9 7-01 
3«0§7f§*0l 
2* 10*00-01 
I • 999 | 3-01 
I •0*710-01 
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§•99*90-02 
3,32092-02 
I • *7 77**02 


VAP1 I ) 
•9 • •*38**0 1 
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•9*22907*02 
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*1 *909«9*02 
•I **7999-02 
*| *999§§*u 2 
-1*39929 -u2 
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•I .00*79-02 
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-3*9*31 0*03 
-2 *§3372*03 
•2*|289|«03 
• l *72*99-y3 
-••0*971-09 
•9* l9«87-U9 


TABLE 2d. Typical pyrotechnic device response data (sample computer printout) 
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ORIGINAL PAGE IS 
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TABLE 2e- Typical meteoroid M-V case positional box score 
(sample computer printout) 


H«s> »*»«P t VCkOc. CROUP | POSITION no, «| 
PR* »j0t7 

Thc RENDON POSITION VECTOR is v fA0t I.OiS 4t«st 
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A VC * ESCAPE VEL* *27972 


TABLE 2f. Typical pyrotechnic device event positional box score 

(sample computer printout) 


PTRQ NO* I POSITION NO* H 
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TABLE 2g 


. Typical recontamination event summary and trajectory 
(sample computer printout) 

HASS feROUf 4 VELOc • 6 R0UP 2 POSITION No. |7 

PR- . | 327 
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)6 

0 
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2 
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TABLE 2g. Typical recontamination event summary and trajectory 

(sample computer printout) (contd) 


STEP 

PATH LENGTH 

X-POSITlON 
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E( r» 
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TABLE 2h. Box score example for all positions of the test meteoroid M-V group 

(sample computer printout) 


•ox score FOR one heTeo*oid of h a ss group 7 and veloo group 

SUhHED OVER 100 POSITIONS 


NO. SAPEHITS 9032.0 
NO. RECONTAH h|TS .*1*30 
NO. NEA* HISSES 3S.I97 
NO. ESCAPES 9 1 • So . 

AVE. ESCAPE VEL. .24S71 


TABLE 2i. Typical final box score and angular escape distribution 
(sample computer printout includes 4 M-V groups; 


•ox score for overall ensemble of me t eonoid mass a no velocity groups 
total no. OF HETEOROIOS in ENSeH&lE DURING MISSION *9 M2. 

NO. SAfChiTS S8 3 * 02 
NO. RECONTAH HITS 3.02b* 

NO. NEAR HISSES IQ3.9§ 

NO. ESCAPES b*8 • B8 


STANDARD OEvUTION FOR NO. OF RtCON T A H i N A T I ON NITS |S 1.7399 
The AVERAGE ESCAPE VELOCITY is .20S0 M/SEC 


The angular distribution for Th* escaped velocity vector relative to th* Direction 

DEG. NO. ESCAPES 

lot 31.1 

20 • 1*70 

3Q* i t • 6 

90* 3S . 1 

SO * 1.17 

* 0 . .*08 

70* 1.97 

• 0* 9|.3 

90» 4*27 

100* 4*99 

llOt 4.99 

120* 4* SO 

130* 30. 9 

Hot *2.3 

1 bo * *1*9 

1 * 0 * 88.1 

170* 79.7 

180 * 89 . * 


RtCON T AHlNATjUN PROBABILITY IS .9&197 
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IV. MODEL EXECUTION 


A. DISCUSSION OF THE INHERENT ASSUMPTIONS 

The total integration of the diverse internal modeling efforts naturally 
requires basic and limiting assumptions. We shall separate these categori- 
cally into those the model unalterably depends upon and those the user may 
arbitrarily define to correspond to a specific mission. 

1 . Fixed Characteristics 

The deterministic meteoroid environment applied to the simulation 
follows precisely the analytical expressions for the fluence reported in 
Ref. 21, a standard and reasonable description. The model used covers the 
interplanetary region of from 1 to 30 A. U. for cometary meteoroids. The 
distribution on the fluence was divided into 10 discrete mass groups for this 
simulation so that a deterministic evaluation would cover the entire range of 
possibilities in 10 passes. This approach avoids the larger number of trial 
meteoroid impacts, necessary for statistical confidence, required during a 
Monte Carlo sequence on this distributed variable. The meteoroid mass 
range is fixed to cover 10" to 10" kilograms because the fluence for 
larger meteoroids fell to below the 1% likelihood range and the smaller 
meteoroids would not disturb the surface enough to eject any dust grains into 
a trajectory. 

The next assumption in this category concerns the spacecraft configu- 
ration. As seen earlier, the spacecraft consists of an element with four 
panels and a separate element. The Z axis is fixed along the Sun vector 
(for maximum photoelectric solar cell usage) allowing computational ease in 
defining shaded areas on the spacecraft. The solar panels are infinitely 
thin-plane sections. The region of the spacecraft considered "sterile" is the 
shaded cone surface and bottom plate of the separate element B (Fig. 15). 

This may be changed in CHKHIT (Appendix C). Any grain, ejected from a 
nonsterile surface, contacting or penetrating these surfaces, is labeled as 
a recontaminating contact. The nonsterile areas, i. e. potential grain ejection 
positions, can be modified in RNDPOS (Appendix C). 
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The interconnecting superstructure of the multiple element body has 
been ignored, for lack of a simplified characterization. 

The spacecraft surface has a residue of particulate contamination that 
is termed individually in space as ejecta. Studies (Ref. 3) have shown that 
the physical dimensions of the dust grains that are typically dislodged are 
distributed predominantly between 10 and 100 |jm, a fixed size range. 

Whenever possible within the coding, proposed trial trajectories would 
first be tested to see if the magnitude of the ejection velocity was greater 
than five times the escape velocity established by the radiation pressure 
force in the direction of the Sun (+ Z axis). Those ejecta having sufficient 
velocity to escape were labeled as such and their trajectory calculation was 
bypassed. 

The grains were considered to be launched along the surface normal 
vector. In the case of the solar panels, an additional 10% of the velocity 
magnitude was randomly added to either the x or y velocity component as 
compensation for the fact that the panels are, in reality, not flat, but struc- 
tured to allow real grain releases along other than normal angles. 

The last fixed characteristic of this model requiring elucidation is the 

particular method applied for calculating the electric field. Surface material 

on this general style of space vehicle consists of predominantly aluminum 

(Al) and silica (SiO-). The thermal blanket covering the element A section 

£ 

was consequently presumed Al; the solar panel illuminated portion was con- 
sidered SiO^, with Al used for the dark side, and SiO^ used to simulate the 
insulating paint coated on the element B exposed surface. Figure 12 shows 
the cross-sectional view of the spacecraft with the assumed component 
materials designated. Other arrangements are possible, but the present 
model can only employ Al and SiO^. 

The electric field in the region of space near any particular surface 
was obtained from the flat-plate solution (Ref. 19) previously discussed with 
the appropriate shaded or illuminated Debye length and surface value from 
the description above. Figure 12 shows, by arrows, the general nature of 
the electric field near the surfaces and in the far field. Three exceptions 
exist: (1) for positions near (compared to the appropriate Debye length) the 

solar panel, the algebraic sum of the surface potential of the illuminated and 
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dark sides is formed and then the field is calculated analogously; (2) for 
positions in the wake of the spacecraft the field is approximated following 
the analysis by Al 1 pert et al. (Ref. 17); and (3) the far-field case is treated 
as an equivalent sphere. 

The discontinuous boundaries unavoidably present in this hybrid 
approximation in a strict mathematical sense are unacceptable. This prob- 
lem comprises the major limitation on the overall validity for the field value. 
Essentially, we consider the grain to be influenced by the characteristics of 
the nearest spacecraft surface component. This assumption proved usable 
upon scrutiny of the initial trial particle trajectories in which the inertial 
characteristics more heavily dominated the movement from step to step, 
whereas the radiation pressure and electrostatic force effects were manifest 
over several intervals of the trajectory. 

2. Mission Peculiar Specifications 

The specific geometric dimensions utilized follow the configuration of 
the multi-element spacecraft and are detailed in Appendix A within the sym- 
bolic computer code element named DATAIN. This element forms the user 
entry point wherein the actual material and physical parameters needed by the 
entire model are defined. The major limitation considered for the surface 
impact submodel is that all the impacted surfaces were treated as 1/8 -inch 
aluminum plate material. The option exists to specify thicker or thinner 
plates of different materials (for instance, the thermal blankets or the solar 
panel silica) via DATAIN, although the choice of the aluminum yielded con- 
servative ejection environments (i. e. , more surface activity). 

The ejecta were considered as silica dust grains with an original sur- 

5 2 

face load areal or density for grains larger or equal to 5;jm of 5. 4 x 10 / m 
(Ref. 22). The integral distribution of areal density has been assumed to have 
an inverse dependence on the square of the particulate size. This dependence 
is set by a DEFINE procedure in RELEAS (Appendix C). Any other decreasing 
function may be substituted, and the normalization is automatic. 

B. OUTPUT DATA PRESENTATION 

As the simulation progressed, evaluations were made, at each meteor- 
oid mass and velocity group level, of the number of impinging meteoroids 
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(denoted as F(I, J) for the (I, J)th mass, velocity group), the number of safe 
hits (S(I, J)) per impacting meteoroid, the number of recontamination contacts 
(R, (I, J)) per meteoroid, the number of ejected grains that escaped (E(I, J)) 
per meteoroid. The mass ranges and velocity ranges corresponding 
to the meteoroid group indices (I, J) are shown in Table 3. 


9 

Fluence = 

1=1 




(mass 

distribution) 


(velocity 

distribution) 


Safe hits 


9 3 

F(I, J) X S(I, J) 

1=1 J=1 



( 22 ) 


(23) 


where S(I, J) is the number of dust grains safely recontacting the spacecraft 
per meteoroid of type IJ as accumulated during simulation. 


Recontamination hits 


9 3 

y y Fd, j) x r(i, j) 

i=i j=i 


(24) 


where R(I, J) is the number of dust grains that recontaminate the spacecraft 
per meteoroid of type IJ. 


Escapes 


9 3 

y y Fd, j) x e(i, j) 

i=i j=i 


(25) 


where E(I, J) is the number of dust grains that escape the spacecraft vicinity 
per meteoroid of type IJ. The overall average velocity of the escaping grains 
is given by: 
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9 3 


Overall average velocity 



A (I, J) X F(I, J) X E(I, J) 


1=1 J=1 


9 3 


(26) 



F(I, J) x E(I, J) 


1=1 J=1 


where A(I, J) is the IJ group average escape velocity magnitude. 

We note that mass-group 10 meteoroids had a fluence of less than 0. 01 
events per mission, and hence, were considered insignificant. 

The summary of the hypothetical full simulation results is shown in 
Table 4. This table has been separated by velocity groups; each group is 
followed by a mission subtotal. The complete enumeration of the model pre- 
dictions is referred to by ensemble grand total on this chart. 

A less conservative interpretation may be extracted from this data by 

assuming only a single nominal velocity for the meteoroid impacts, for 

-4 

instance, velocity group 1 with the range 1.2 - 2.0 x 10 m/s, and linearly 
interpolating (multiply by 3 to simulate the number of total trials) the results. 
Also considered here was a pyrotechnic event simulation (see Appendix B). 
Comparison of these three approaches appears in Table 5. 

The data from Table 4 have been graphically plotted in a perspective 
3 -dimensional manner as seen in Figs. 16-20. The most immediately 
noticeable characteristic is the strong inherent velocity dependence. These 
plots all have linear Z-axes and the peak value indicated. The anomalous 
hump at the low end of the fluence surface exists because the first mass 
range increment is larger than the remaining intervals (see Table 3). It is 
apparent that the small meteoroids (10 ^ -3 X 10 - ^ kg) heavily dominate 
this distribution. The next figure represents the distribution of grains that 
were recaptured safely by the nonsterile portion of the spacecraft. It is 
interesting to note that the three surfaces for safe hits, recontaminating 
hits, and escapes, respectively, contain peaks within their range of validity. 
This means the competing aspects throughout the complete simulation appar- 
ently did span the range of actual physical interaction within the limits of 
the assumptions. 
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The re contamination data are shown in Fig. 18. The nature of the 
double hump is elusive at present. The spacecraft may have areas that 
could be considered as M hot spots. !? The random selection of the metedroid 
impact position considered 100 choices for each mass-velocity group. The 
bulk of the recontaminating events were from either the top of the solar 
panels or the element A thermal blanket material. (This is an area that 
would most likely yield lower predictions for the recontamination, should the 
element A surface be assigned a more accurate thickness and other mechani- 
cal property values. ) Perhaps the solar panel and element A surfaces have 
different "pass bands" for the recontamination events. The surface posi- 
tions allowed for meteoroid impact were on the illuminated region of ele- 
ment A, both sides of the solar panels, and the illuminated cone of element B. 
Several surfaces were positively charged, requiring the ejecta to initially 
have a positive charge. The sterile zone is always shaded and therefore has 
about -60 eV in electric potential. When the ejecta were forced toward ele- 
ment B by the solar radiation pressure, there were perhaps a few ideal pos- 
sible trajectories allowing each of these areas to bring the positive ejecta 
close enough to element B to be attracted. 

Figure 20, for the average escape velocity, indicates that meteoroids 
-8 

of up to 3 X 10 kg and for all velocity ranges cause grains to escape with 
about the same average velocity, and not until the large meteoroids impact 
(groups 8 and 9) does the typical value of about 0. 2 m / s drastically change. 
The relative importance of these escape velocities was studied by weighting 
the average velocity of any particular M-V group by the number of grains 
escaping, and then normalizing by the total number of escapes. That is 


Weighted velocity magnitude (I, J) = 


A (I, J) * F(I, J) * E(I, J) 
9 3 

7 7 F(i > j) * e(i * j) 

1=1 J=1 


(27) 


Figure 21 indicates this functional surface with the high mass -velocity 
groups dominating the spectrum of escaping grains velocities. 
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During the execution, as mentioned in Section III. C, the particles 
escaping were classified and summed into 18 angular bins corresponding to 
0-180 degrees in 10-degree increments. Figure 22 contains the polar plot 
of the resultant angular distributioh obtained during the simulation. The 
spikes in the Sun direction may be explained by the predominant number of 
grains that exceeded the initial escape velocity criteria from: 1) the top 

plate and upper solar panels, 2) the illuminated portion of element B, 
and 3) the element A cone surface. The bulk of the grains were forced 
outward from the Sun by the solar radiation pressure. This polar plot lends 
further confidence to the overall model and indicates that the electric field 
approximation was realistically bounded by the momentum properties (the 
spikes 1, 2, and 3 on Fig. 22) and the solar radiation pressure force (area 4 
on Fig. 22). 

C. OPERATIONAL ECONOMICS 

This model was executed on a UNIVAC 1108 computer system in a 
sequence of runs. The compiled processor occupies approximately 16K of 
core and the Monte Carlo/ deterministic structure is inherently CPU bound. 

Each of the M-V groups and 4 pyrotechnic firings were evaluated 

totaling 18. 5 hours of CPU time. The present computer rates at JPL 

5 

amounted to approximately $4. 7K for a grand total of 2. 5 x 10 histories 
(some runs were made during Prime Time at increased expense). This 
gives a rough factor of about 2^/history. 
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Table 3. Discrete mass and velocity group 
index -magnitude correspondence 


Mas s 
group I 

Mass 

range, kg 

Velocity- 
group J 

Velocity range, m/s 

1 

10- 15 

-14 

-3X10 

i 

1.2 - 2. 0 X 10 4 

2 

3X10- 14 

- 3 X 10" 13 

2 

2. 0 - 2. 8 X 10 4 

3 

3X10- 13 

- 3 X 10" 12 

3 

2. 8 - 3. 6 x 10 4 

4 

3 X 10" 12 

- 3 X 10 -11 



5 

3 X 10" 11 

- 3 X 10 -10 



6 

3 X10" 10 

- 3 X 10‘ 9 



7 

_9 

3 X10 7 

- 3 X 10" 8 



8 

3 X 10 -8 

- 3 X 10" 7 



9 

3 X 10 -7 

- 3 X 10‘ 6 
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Table 4. Output data distribution from the simulation of the hypothetical spacecraft mission 
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Table 5. Summary of results of recontamination analysis 
for a hypothetical spaceflight 
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Fig. 16. Mission micrometeoroid fluence distribution 
over meteoroid mass and velocity groups (after NASA 
Ref. 19) 
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Fig. 17. Mission safe-hit distribution over meteoroid mass and 
velocity groups 


JPL Technical Memorandum 33-737 








Fig. 18. Mission grain recontamination distribution over meteoroid 
mass and velocity groups 
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Fig. 19. Mission grain escape distribution over meteoroid mass 
and velocity groups 
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PEAK AT (9,3) 
OF 4,94 m/s 



Fig. 20. Mission average escape velocity distribution over meteoroid 
mass and velocity groups 
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Fig, 22, Escaped grain 
angular distribution 
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V. CONCLUSIONS 


The major objective of this research was to study the effects of typical 
mission environments on the distribution of particles on spacecraft surfaces. 
Specifically, the various migration mechanisms that result in particulate 
redistribution were investigated and quantified. The ultimate goal was to 
develop a methodology and quantitative analytical tool for the evaluation of 
the recontamination hazard for various planetary missions and mission 
strategies. 

The results obtained for a simulated hypothetical mission of 300 days 
are given in Table 5. These indicate that within the limits of the stated 
assumptions, a total of 28, 279 meteoroids impacted the craft causing 148, 026 
grains to be ejected (less than 1% of the total load) with 126, 778 escaping, 

21, 215 recontacting the spacecraft on safe areas, and 33 recontaminating 
element B sterile zone. The pyrotechnic events added another 35, 000 
escapes, 849 safe hits, and 1 recontaminating grain. Assuming a Poisson 
distribution for the statistics, we see that the expected number of particulate 
recontamination events is 33 with a standard deviation cr ( 33) = SQRT(33) = 5.7 
and that the probability of one or more recontaminating events occurring is 

Prob (recont. ) = 1. - exp(-33) = . 9999 . ♦ . (28) 

The numbers predicted scale directly with mission elapsed time and particu- 
late areal density. 

The general structure of the computer software (see Appendix C) 
allows the ready application (with few modifications) to studies of varying 
natures. For instance, the sterile zone may be redefined easily to be the 
cylindrical area of space observed by the Canopus sensor in the near field 
to count the predicted number and velocities of the grains traversing its 
view; particular meteoroid classes may be studied; various impact locations 
may be designated; spacecraft geometries that may be interpreted as degen- 
erate cones are applicable; the ejected particles captured by planetary bodies 
during flyby may be studied; etc. 
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Each particular spaceflight mission simulation requires specific 
parameter definition in order to identify material constituents, interplane- 
tary setting, etc. The following computer listing. Table A-l, shows the 
symbolic code used for the hypothetical mission discussed earlier. 

The common block section provides the interconnection map which 
transfers data throughout the general computer program. The parameters 
have been arranged into groups corresponding to their general use. The 
first grouping is for the meteoroid impact model and grain release routines. 
The variable ACCMIN is the threshold level for surface acceleration below 
which dust grains are not released. Array ADRHO contains the material 
density of the ejected particle for aluminum, graphite, and silica. AKM and 
SIG are parameters that were obtained from a best fit of the particle adhesion 
model to an experiment (Eq. 4). APVEL is the array that holds the three 
meteoroid velocity group magnitudes. E is Young's modulus for aluminum. 

H is the presupposed thickness of the target surface. MAT is the variable 
that designates the material assumed for the ejected dust grain. PO is the 
loading function (detailed analysis appears in Ref. 3). Finally, RHO speci- 
fies the target density (aluminum). 

The next grouping concerns the parameters necessary for electric field 
and potential determination. The AMAT array determines the material 
(either aluminum (Al) or silica (SiC> 2 )) for the spacecraft sections as indi- 
cated. The arrays A and B contain material dependent coefficients obtained 
from an empirical equation for the secondary-emission coefficient of the 
form: 


A (I) ' B (I) e *o 

kT - e<b 
e 


(Al) 


where 1=1, 2, 3 for aluminum, carbon, and silica respectively, kT e is the 
Boltzmann constant times electron temperature and e<\> Q is the current poten- 
tial (both in eV). (Note: Actually A and B also depend weakly on kT e , so 
values are valid only for kT e near 20eV. ) Array EL is the Hollerith dic- 
tionary reference for the AMAT array. The other electrical parameters 
are self-explanatory. 
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The next grouping covers the basic spacecraft parameter identification 
scheme. The origin is one Debye length from the bottom plate of the model. 
The parameters SP1-SP4 and SP01-SP03 are measures of distances in' the X 
axis direction, whereas D1-D4 indicate the Z axis dimensions of the cone 
altitudes and inter-element spacing. Following the geometric inputs are 
arrays that contain solar spectrum energies and yield values as noted. 

Finally, the miscellaneous variables appear and are mostly self- 
explanatory. K1 must be fixed at 10. The entries for N1 and N2 end up 
being default values because they are read in again via input case card. 

SRMIN and SRMAX specify the heliocentric limits for the mission, whereupon 
the fluence evaluations are computed using the geometric mean of these end 
points. 
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TABLE A-1. ELEMENT DATAIN SYMBOLIC LISTING 


c 

c 


c 

c 

c 

c 

c 

c 

c 


c 

c 


r 


c 

c 


k 

c 

c 


c 

c 


c 

C- 

c 

c 

c 


c 

c 

c 

c 


BLOCK DATA 8 INFUT OATA FOR THE SUBROUTINES. 

COMMON BLOKCl CONNECTS DATAIN. MAIN AND RELEAS. 

COMMON /BLOKCl/ AKM » DRHO »K1»’R9VEL(35»10) ,CR(1G)«SIG« OOIA(IO) 

It ANP ART.ANNORM .PRBDIA(IC) 

COMMON BLOK C 5 CONNEC’S OATAIN, MAIN, RELEAS AND Y ANSI 
COMMON /BLCKCS/ KO »PKASS ,PRH C»RS ( 15 ) »FV AP ( 3E J , FPK ACC ( 35 > t RHO 
COMMON BLOK OB CONNECTS DATAIN, MAIN ANO YANG1. 

COMMON /BL0KC6/ E»H.NU»PO»ISKIP 

COMMON BLOKOB CONNECTS CONTAM, CAT AIK 1 , MAIN AND PARPOT. 

COMMON /BLOKC9/ DIM , PATHMX ,SR » VELMA C 

COMMON BLOK 17 CONNECTS DATAIN, MAIN ANO YANG1. 

COMMON /BL0K17/ PVEL.ACCMIN 

COMMON BLOK 18 IN DATAIN»EFIELO»ESURF*MAIN»RNCPOS»SCPLOT» SHADE* ANO CONTAM 
COMMON /BLOK 18/ SP1»5P2,SP3»SP4 

COMMON BLCK19 CONNECTS CATAIN, EFIELO, ESU RF .MAIN, RNOPOS ANO SCPLOT. 
COMMON /BLOK 19/ SPC1.SPC2.SF C3 

COMMON BLCK24 CONNECTS E e IELD , ESURF MAIN AND DATAIN. 

COMMON /BL0K24/ALAM8(2) * ANPHCT (Cl, AKFHOT <2 ), EFEL ( 2 J , EFED ( 2 ) ,DEB . 

1AM ATI 9 J ,EH2) fPHIEL «2) »PHIEC< 21 . A LAMA V. PHIAVE 
COMMON BL0K2E CONNECTS. DATAIN, MAIN AND PARPCT. 

COMMON /BL0K26/ MAT.AVR 

COMMON BLCK2 S CONNECTS DATAIf. AND PARPCT. 

COMMON /BLOK 2 8/ A (3) ,B ( 3 > » MATP OS ( 3) ,S OLS PC( 13 5 . 2 ) , Y IELO <15, 3 ) , Y IEL 
1DPC16»:> 

COMMON BLOK 29 CONNECT OAT AIN .EFIEL3, ESU RF . MAIN ANC PARPOT. 

COMMON /EL0.K29/ AKTE.AKTP.ANE 

COMMON BLOK3C CONNECTS *DA T AIN AND MAIN. 

COMMON /BLCK3C/ADRHO (3> »D1 »D 2, 03 , C4 »Nl » N2 , SRMAX.SRMIN 
1. EPSILN, ELECT ,APVEL C3) .TMISS , KPLOT 

INPUT DATA FDR METEOROID IMPACT AND GRAIN RELEASE ROUTINES. 

DATA ACCMTN/D .0 31/ a MINIMUM ACCELERATION OF INTEREST. 

DATA ADRH0/?7rC.»225C.»2RC0./ -3 ALUM., GRAPHITE AND SILICA. 

DATA AKM/-.882/ SIG/.400/ 8 STATISTICAL QUANTITIES NEEOED FOR RELEAS. 

DATA APVEL/l.EE4,2.4E4,3.2r4/ 8 METEOROIC VELOCITY VALUES. 

DATA E/6.397E10/ 3 YOUNG’S MODULUS. 

CATA H/3.175E-3/ 3 TARGET THICKNESS (METERS.) 

MAT -- NUMBER DEFINING MATERIAL OF PARTICLE. 1 = ALUMINUM, 

2 = GRAPHITE, 3 = SILICA 
OATA MAT/ 3/ 3 CHOSEN FOR SILICA. 

DATA NU/.25/ 8 POISSON’S PATIO. 

OATA P0/1.74E9/ 8 LOADING FUNCTION I N/M*M. ) 

DATA RH0/2.G83E3/ S TARGET DENSITY. (ALUMINUM.) 

GENERAL INPUT CATA FCP THE ELECTRIC FIELD AND SURFACE POTENTIAL 
ANO PARTICLE POTENTIAL ROUTINES. 

AMAT SPECIFIES THE MATERIAL OF THE SPACECRAFT SECTIONS. 

DATA AMAT/’AL •» *AL •» ’SIC2* , *AL* ,*AL*» *AL’, *SIC2»» *SI02*» ’SI02’/ 

DATA A/7. 67, 2 . 76. 10. 17 / 3 WEAK DEPENDENCE ON AKTE. 

DATA B/ .430, .122, .382/ 8 WEAK CEPENDENCE ON AKTE. 

DATA EL/*AL*.’ST02’/ 3 FOR ALUMINUM AND SILICA. 

DATA AKTE/20 ./ 3 BCLTZMAN CONSTANT TIMES THERMAL ELECTRON TEMP. 

DATA AKTP/10. 7/ 3 BOLTZMAN CONSTANT TIMES THERMAL PROTON TEMP. 

DATA AVP/4.2EEE/ 3 THE PROT CN DIRECTED VELOCITY. 

DATA ANE/5.E6/ 3 PROTON NUMBER DENSITY. 

DATA AKPHOT/3 . .2 •/ 3 PHOTOELECT RON TEMP (EV) 

DATA ANPHCT/1 .E9.1.E8/ 8 PHOTON NUMBER DENSITY. ( MKS ) 

ANPHOT/ANC .GT . AKTE/AKPHOT .GT. l.C MUST OBTATN 


SPACECRAFT DIMENSIONS 
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THE ARRAY SOLSPC CONTAINS VALUES OF ENERGY IN ELECTRON VOLTS < TRANS- 
LATED FROM WAVELENGTH OF THE SOLAR SPECTRUM) AND THE SOLAR 
IRR ADIANCE AS AN INTEGRAL FR CM M AV ELENGTHzC . C CR E=INFINITY. 


DATA 

(SOLSPC(I.l) , 

1=1. 62) 

/ 




1 

.00 

124. .01240. .01550. .02067. .02480. .03100 

2 

.0 3543 » .04133. 

•04S60 * 

•C62C0, 

.06526, 

.06889, 

.07294, 

3 

.07750. .08267. 

.08857, 

. 08538, 

.10333, 

.11273, 

.12400, 

4 

.13778. .15500. 

•17714, 

•2C6 67, 

.24ecc, 

.25306, 

.25833, 

S 

.26383. .26957. 

.27556, 

.28182, 

.28337, 

.29524, 

.30244, 

8 

•31CG0 . .31795. 

•32632, 

.33514, 

•34444 , 

.35429, 

.36471, 

7 

.37575. .38750. 

.40000, 

.41333, 

.42759, 

•44286, 

.45926, 

8 

.47692. .49600. 

• 51867 f 

.53913, 

.56364, 

.59048, 

•6200C , 

8 

.65263. .68889. 

.72841, 

.77500, 

.32667, 

.88571, 

.95385 / 


DATA (SOLSPCd.l) . 1=63. HE) / 

5 1.03333. 1.12727. 1.24003, 1.30526. 1.37778. 1.45882. 

1 1.550CC, 1.65333. 1.6756E. 1.69863. 1.72222. 1.74648. 

2 1. 77143. 1.79710. 1.82353. 1.85075. 1.87879. 1.90769. 

3 1.9375C. 1,96825. 2.0C0CC. 2.C3279. 2.G4959. 2.C6667. 

4 2,09403. 2.10170. 2.11965* 2.13793. 2.15652. 2.17544. 
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TABLE A-1. ELEMENT DATAIN SYMBOLIC LISTING 


5 

2.19469* 

2.21429* 

2.23422* 

2 »25<tf5» 

2.27523* 

2.29630* 

6 

2.31776* 

2.33962# 

2.36190* 

2.33162» 

2 .40777 * 

2.43137* 

7 

2.45545# 

2.48000* 

2.50505* 

2.F3rEl» 

2.55E7C* 

2.56333* 

8 

2.61053* 

2.63830* 

2.66667* 

2.63565. 

2.72527* 

2.75556 / 

DATA CSOLSPCf I»l> » I 

=117* 170) / 



9 

2.78652# 

2.81818* 

2.35057* 

2.88372* 

2.31765* 

2.35238# 

5 

2.98795* 

3.02439* 

3.C6172* 

3.1CCCC* 

3.13924* 

3 • 17949 * 

1 

3.22C78* 

3.26316* 

3.30667* 

3.35135* 

3.33726* 

3 .44444# 

2 

3.49296* 

3.54286* 

3.5942C* 

3.64706* 

3.7C145* 

3.75758# 

3 

3.81538* 

3.37500* 

3.93651* 

4.00000* 

4.06557* 

4.13333* 

4 

4.20339* 

4.27586* 

4.35088* 

4 .42857* 

4.509C9* 

4.55259# 

5 

4.67925# 

4.76923* 

4.36275* 

4.96000* 

5.06122* 

5.16667* 

6 

5.2766C# 

5.39130* 

5.51111* 

5 .63636* 

5.90476# 

6 .20000 * 

7 

6.52632* 

6.88889* 

7.29412* 

7.75000* 

8.26667* 

9.85714 / 


DATA CSOLSPC(I#l) t 1=171 * 1851 / 

8 10* 33333* 10*7 8261* 11.27273* 11*80952* 12.40000* 13*05263* 

9 13 *77778 # 14.58824* 15.5CCC0* 16.53333* 17 .71425* 19.C7692* 

5 20.66667* 22.54545* 24.83000 / 

DATA 1 SOLSP C ( 1*2 1 * 1=1* 361 / 125. 3C# 135.29865* 135.29855* 

.29735* 135.29671* 135.29556* 
,28596 * 135.28376* 135 .28101* 
.26801* 135.26091* 135.25056* 
>18356 » 135.13456* 135 .05506* 


1 

2 

3 

4 

5 


135.29829* 

135.29328* 

135.27756* 

135.23606* 

134.91806* 


135.29801* 
135.28776* 
135.27321# 
135 .21506 * 
134.63906* 


135 

135 

135 

135 

134 


6 

134.403* 

1 34 .3*1 1 » 

134.272* 

134 ,199. 

134.116* 

134.025* 

DATA 

« S0LSPCII*2). 1=37* 90) 

/ 



7 

133.819* 

133. 7C2# 

133.572* 

133.433* 

133.277* 

133.098* 

3 

132.889# 

132.646* 

132.361* 

132.031* 

131.661# 

131.251* 

9 

130.796# 

130.286* 

129.696* 

129.036* 

128.301* 

127.456* 

5 

126.491* 

125.346* 

123.921* 

122.116* 

119.386* 

117.231* 

1 

114.116* 

110.456* 

106.056* 

ICC. 676* 

93.9655* 

90 .0334* 

2 

85.7234* 

81.0309* 

75.7934* 

69.9384* 

68.6909* 

67.4159* 

3 

66.1139* 

64.7849* 

83.4284* 

62.0429* 

60.6284* 

52.1869* 

4 

57.7159* 

56.2174* 

54.6893* 

53.1329* 

51.5469* 

49.9284* 

5 

49.1079* 

48.2797* 

47.4421* 

46.5572* 

45 .7442* 

44 .8874 / 

DATA 

( SOLSPC 11*21* I 

- 31 * 144) 

/ 



6 

44 .0289* 

43.1712* 

42.3165* 

41.4689* 

40 .6132* 

39 .7519# 

7 

33.8822# 

37,3979* 

37.P377* 

35.1327* 

35.2592* 

34.3379* 

8 

33.4214* 

32.4927. 

31.5422* 

3C.5767* 

29.6012* 

28 .6237* 

3 

27.6422* 

26.6237# 

25.6002 * 

24.58C9* 

23.5607* 

22.5322* 

5 

21 .5014, 

20.4857* 

19.5C27# 

18.57T7* 

17 .7C2 4 * 

16 .8769* 

1 

16.0439* 

15.1333* 

14,3037* 

13.4224* 

12.5737* 

11.8054* 

2 

ii .i5rs* 

10.5792* 

10,0302* 

3.47566* 

9.90641* 

8.32191* 

3 

7.74365* 

7.19366* 

6.65591 * 

6.11191* 

5.57141* 

5.03566# 

4 

4.49691* 

3.96191* 

3.45341* 

3 .00 2 IF * 

2.60366* 

2.24041 / 

DATA 

(SOLSPCfI.21* I 

= 145. 1351 

/ 



5 

1. SI 74 1 * 

1 ,63816* 

1.36266* 

1 .03716* 

.857910* 

.763660* 

6 

.657160* 

.543160# 

.443913 . 

.365160* 

.306660* 

.263060* 

7 

.227385* 

.193560* 

• 16298 C » 

.131485* 

.096585* 

.067985# 

3 

• C27 78 5 • 

• C 10 935 * 

.004232 . 

•C023CC* 

.CC135C0 

* .000930# 


9 

5 

1 


• GDO 7 8 • 

. CQO 5372 1 
.0005010* 


.CCC73* • COC 5999 3 
• QCC 5292 * 

.0 DC4 984 * 


599* 134.556* 134 


510* 134.459* 


.0005215* 

.CCC4929* 


.00057200* .CCC5S57C* .0005526* 
•0C05121. .0005059* .0005027* 
• CCC4896 * .0004863 ✓ 


C 

c 

c 

c 


THE ARRAY YIELD CONTAINS TH6 VALUES 
range FOR THREE MATERIALS. 


FOR THE YIELD IN EACH ENERGY 


DATA YIELD 
6.3E-3* 
2.7E-7* 
4.7E-3* 
0. * 0.9 
8.5E-2# 


/l.CE-7# 3.1E-6* 3.GE-6* 


2 

3 

4 

5 

6 


3.2E-2# 
1.1E-5* 
1. 8E-2* 
l.CE-4* 
.1 * • 1 * 


7.4E-2* 
3.2E-5* 
2.5E-2* 
1.9E-3* 
.063 / 


.15* .2* .23* . 
3.0E-5* 3.1E-4* 
2.5E-2* 2.5E-2* 
l.OE-2# 2.302* 


3.1E-5* 1.6E-4* 3.6E-4* 

25* .27* .12* C.* C.« 
7.4E-4* 1.3E-3* 

2.5E-2* C.C* C.* C.* 
4.0E-2* 6.8E-2* 
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TABLE A-1. ELEMENT DAJAIN SYMBOLIC LISTING 


THE ARRAY YIELDP CONTAINS: FIRS T t SEVERAL VALUES OF ENERGIES 

• ELECTRON-VOLTS) DIVIDING THE YIELD « ELECT RONS /PHOT ON ) OATA INTO 
RANGES. AND SECOND THE AVERAGE ENERGY VALUE FOR EACH RANGE. 

DATA YIELDP / 4.0. 4.25. 4.75. 5.5. 6.5. 7.5. 8.5. 9.5. 

2 10.5. 11.5. 12.5. 13.5. 14.5. 15.5. 22.5. 27.5. 4.12. 4.5. 

3 5.. 6*. 7*. 8 • . 9*. ID*. 11.. 12.9 13.9 14.. 15. * 19.. 

4 25 .. 0. / 

MISC. INPUT FOR THE PAIN DRIVER. 

DATA ANPART/E.4E05/ 8 CONTAM AREAL DENSITY. GT. 5 MICR0NSCM**-2 ) 
DATA ELECT/t .6E-19/ 8 ELECTRONIC CHARGE <MKS.) 

DATA EPSILN/e.e5E-12/ 8 FREE SPACE PERMITTIVITY. 

DATA KPLOT/O/ i O-NO PLOT. 1-PLOT THE SPACECRAFT. 

DATA K1/1C/ S NUMBER CF EJECTA DIAMETER BINS. 

DATA N1/2/N2/4/ 8 OUTER AND INNER LOOP INDICIES RESP. 

DATA P AT HMX/2C ./ 8 TWENTY METERS FCR WANT CF A VALUE. 

DATA PRHO/500. / 8 METEOROID DENSITY. I MKS I 

DATA SRMAX/1 .5/SRMIN/l •/ 8 HELIOCENTRIC RAOII LIMITS. 

DATA TM ISS/300 ./ 8 MISSION DURATION - DAYS. 

END 
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PYROTECHNIC EVENT MODEL 
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The pyrotechnic model used for the simulation is comprised simply of 
two parts: a plausible value for the surface acceleration related to the particular 
pyrotechnic device of interest (that decreases exponentially with range - see 
example in Table 2), and a means of localizing the pyrotechnic device on the 
geometric model. The choice of the surface acceleration, of course, depends on 
the pyrotechnic device and will not be discussed further. The means of localiz- 
ing the device is accomplished by using the standard random position selec- 
tion subroutines (RNDPOS) modified as follows. In Fig. 15, the names of the 
areas of the model that are considered contaminated are: ATP, the area of 
the element A top plate; ATS, the area of the element A cone surface; ASP, 
the solar panel area (both top and bottom); and All, the nonsterile cone sur- 
face area of element B. A linear mapping relates scalar values (normalized 
by the total ATP + ATS + ASP + All) to geometric areas on the spacecraft. 

The input of two limiting scalar quantities (e. g. , the pyrotechnic data card 
on Table 2) establishes the definition of a specific subarea of spacecraft for 
the pyrotechnic event grain ejection study. 
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The FORTRAN symbolic codes describing the main driver and the 
subroutines are now listed. The brief general nature of each element is 
shown in Table C-l. 


Table C-l. Brief functional description of the driver programs 

and the subroutines 


Program 

Names 

Functional Description 

CHKHIT 

Cone region surface recontact determination. 

CONTAM 

Particle transport control coding (no printout). 

CONTM2 

Same as CONTAM with the trajectory, etc. printed. 

E FIELD 

Electric field of the model spacecraft. 

ESURF 

Spacecraft surface potential calculation. 

FULBPS 

Driver used when surface effect and meteoroid model 
data input on cards. 

FULDET 

Driver that contains the meteoroid model/impact model. 

PARPOT 

Potential of the ejected particle. 

PYROM 

Pyro event driver. 

RANDNO 

Random number generator. 

RELEAS 

Particle release effect model. 

RNDPOS 

Random ejection location code. 

SCPLOT 

Plot profile of spacecraft. 

SHADE 

Determine whether ejecta are shaded. 

SPHIT 

Determine if a solar panel was hit. 

THINPL 

Thin -plate meteoroid impact model. 

TRAJEC 

Incremental particle trajectory code. 

TRANSL 

X-Y Symmetry maintenance code. 

VNORML 

Calculate initial normal vector for ejecta launch. 

YANG1 

Driver for thin -plate impact code. 
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SUBROUTINE CHKHIT 


SUBROUTINE CHKHIT 

C •••• NASA JPL •**• 9/12/74 *•** O.EDGARS* BIONETICS •••# 

♦ • * * 

THIS ROUTINE CHECKS WHETFEP t DEVIOUSLY FREE FARTI CLE (PO S VEC J COMES 
IN CONTACINEWPOS) WITH THE CONICAL SPACECRAFT AREAS. THE SC MOTEL 
CATA REQUIRED APE ACCESSED VIA COMMON BLOCK AND ARE DEFINED AS 
ENCOUNTERED. 

A SEPERATE SUBROUTINE MCNITCFS THE SOLAR PANEL ZONE AND IS CALLED 
WHEN THE PARTICLE APPROACHES THEIR PROXIMITY. 

definitions: 

111 XYZlJI It THE NEWIPRESENTJ POSITION OF THE PARTICLE. 

121 IS IS THE INDEX CENCTING T FACEC RAFT SECTOR OF INTEREST. 

IT) IM1 IS THE SAFE SURFACE HIT INOEX. (OUTPUT! 

IA) r-DID NOT CONTACT SAFE SURFACE 
191 1-5IO CONTACT SAFE SURFACE 

(Ml IH2 IS THE INDEX DENOTING 

CA) D - DID NOT CONTACT STERILE SURFACE 

IB) 1 - DID CONTACT STERILE SURFACE. 

• * •• 


C 


c 

c 

c 


REAL NEWPCS 

COMMON BLOK 02 CONNECTS CHKHIT. MAIN. RNCPOS AND SHAOE. 

COMMON /BL0KC2/ ALI4I 

COMMON BLOK 03 IN CHK HIT » CONTA M, MA IN.S HADE .S PHIT »T RA JEC AND TRANSL. 
COMMON /BLCKC9/ NEWP0SI3) 

COMMON BLOK 13 CONNECTS CHKHIT. CONTAM. MAIN AND SPHIT. 

COMMON /BL0K13/ IHITI2) 

COMMON 3L0K15 IN CHK HIT » CONTA M. EFIELO . MA IN. RNDPOS .SHADE. VNORML. 
CCMMCN /BLGK15/ ISECTR 

COMMON 3L0K2C IN CHKHIT* EFIELD. MAIN. RNDPOS . SHADE AND VNORML. 
COMMON /BL0K2C/ TANTH1.T ANTH2.TANTH7.T ANTH4 

COMMON 3L0K21 IN CHKHIT* CONTA M, EFIELD . MA IN. RNDPOS .SHADE. SPHIT. SCPL 
COMMON /BLOK 21/ Z1.Z2.Z3.Z4.Z5 

EQUIVALENCE I XYZ (1) .X EWPOS 111 I . I ISECTR. IS) .(IHITIII.IHII.I IHITI 2 I. 
1IH2I 

DIMENSION XYZ 1 3 I 
GO TO! 100.2CC. 300 » IS 

100 IFIXYZI3I.GT.Z5I GO TO 15C 3 IN SECTOR 1 

ZSC= SQRTIXVZ(1I*XYZ|1I ♦XYZ (21 *XY? (2) I »T ANTH3 ♦ AL(3) 

IFIZSC.LT. XYZI3II GO TO ISO 
12S IH1=1 
IH 2=0 

RETURN 3 HIT SAFE SURFACE. 

ISO IH1=C 
IH2=C 

RETURN 3 NO HITS REGISTERED. 

20C IF(XYZI3) .GT.Z3I GO TO 150 3 IN SECTOR 2 

ZSC= SQRT(XYZ(1I*XYZ(1I*XYZ(2I»XYZ(2> )*T ANTH2 +AL(2) 

IFIZSC.LT .XVZ(3I) GO TO 150 
GO TO 125 

3CC IFIXYZI3I .LT .Zl) GO TO 150 3 IN SECTOR 3 

ZSC= SQRTI XYZ II I • XYZ (ll*XYZIZI»XTZ(2) )*T ANT HI ♦ALII) 

IFIZSC .G T .XYZI3I I GO TO 150 
IH1=C 
T H2 r 1 

RETURN a HIT STERILE SURFACE. 


END 


ORjgj&AL PAGE B 
OF POOR QUAI*G3[ 
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SUBROUTINE CCNTAK 


C NASA JPL 9/12/74 D. EDGARS* BIONETICS •••• 

c •***• • 

C THIS SUBROUTINE COORDINATES THE TRAJECTORY ANALYSIS FOR THE EJECTED 

C GRAIN UNTIL THE GRAIN ACHIEVES ONE OF THE r OLLOM IN6 CONDITIONS: 

L. 

C 1-THE GRAIN RECONTACTS AN A CCEPTA3LE PORTION OF THE SPACECRAFT 

C WHEREUPON THE VARIABLE IHTll) IS SET TO 1 - SAFE HIT. 

£ 

C 2 -THE GRAIN TRAJECTORY SURFASSES THE PRESCRIBED MAXIMUM PATHLENGTH 

C OF PATHMX (20 METERS* - ESCAPED. 

r 

v 

C 3-THE NUMBER OF TRAJECTORY INCREMENTS EXCEEDS 500 STEPS - ESCAPE. 

C 

C 4-THE GRAIN HAS RECONTACTED AN AREA OF THE SPACECRAFT CONSIDERED 

C STERILE WHEREUPON IHIT { 2 ) IS SET TO 1 - RECONTAM Hit. 

C 

c curing The EXECUTION The particle is KINEMATICALLY PROJECTED FROM 

C AN INITIAL POSITION THRU INCREMENTED POSITIONS USING THE PARTICLE 

C CHARGE* THE ELECTRIC AND INERTIAL FORCES EXPERIENCED BY THE 

C PARTICLE AT THE INITIAL LOCATION OF EACH INCREMENT. THE NEW POSI- 

C TION IS SYMMETRICALLY MAINTAINED IN THE FIRST OCTANT AND MONITORED 

C FOR THE HIT CONDITIONS MENTIONED. 

C 

C THE PARAMETERS ARE 


c 

PCS OLD ( 3 ) 

- POSITION VECTOR 1INITIAL CONDITIONS). 


INPUT 

c 

VELOLDt 3) 

- V FL 0 C I T Y VECTOT 1 INITIAL CONDITIONS). 


INPUT 

r 

GP AS 

- THE PAPTICLF ( EJECTED MASS) 


INPUT 

C 

VMAG 

- THE VELOCITY MAGNITUDE 


INPUT 

c 

ISC 

- DEFINES THE SECTOR OF CONCERN. 


INPUT 

c 

TP 

- DEFINES CONDITION WHERE t 


INTERNAL/OUT 

c 


1 = PARTICLE UNDER A SOLAR PANEL. 

c 


0 = PARTICLE NOT UNDER A 

SOLAR PANEL. 

c 

ISH 

- DEFINES CONDITION WHERE: 


INTERNAL/CUT 

c 


1 r PARTICLE SHADED FROM 

SUN ( OUTSIDE HAKE) 

c 


C = PARTICLE NOT SHADED 

FROM 

SUN. 

c 


-1 = ' PARTICLE SHADED FROM 

SUN 

(INSIDE WAKE) 

c 


• • * •* * 


*•** 


REAL NEWPOS* NCWVEL 

DIMENSION PCSCLD«3)» POSNEW (3 ) » VELCLD(3>» VELNEW € 3 X 
C COMMON BLOK 04 CONNECTS CONTAM* MAIN* PARPOT* SPHIT AND TRAJEC. 

COMMON /BL0KC4/ COEFPP » CON ST A* fRTCH£*PHIO 
C COMMON BLCK03 CONNECTS CONTAM. DATAIN* MAIN ANO PARPOT. 

COMMON /BLOKCS/ DIM . F AT H MX *S fi . VELM A C 

C COMMON BLCK09 IN CHK HIT » CO NTA M, MA IN, S HADE *S PHIT »T RA JEC ANO TRANSL. 

COMMON /BL0KD9/ NEWPCS(3) 

C COMMON 3L0KIC CONNECTS CON T AM . M AI N*S PHIT . T R A JEC ANO TRANSL. 

COMMON /BL0K12/ NEWVEl(3>» DELTAT 

C COMMON 3LCKU IN CON T AM » EFIELD . MA IN. R NOPOS *S PHIT *T RA JEC. V NORHL. 

COMMON /BLOKU/ PCSVEC(3» 

C COMMON BL OK 12 CONNECTS CONTAM, MAIN ANO TRAJEC. 

CCKMON /BLOK 1? / CM ASS *VELVEC(I> »ILCST 
C COMMON BLOK 1 1 CONNECTS CHKHIT , CONTAM, MAIN ANO SPHIT. 

COMMON /BL0K1?/ IHIT ( 2 ) 

C COMMON 9LCK14 CONNECTS CONTAM, EFIELC, MAIN, SHAOE ANO TRAJEC. 

COMMON /OL OK 14 / Ir AN EL » ISH AC E 

C COMMON BL0K1B IN CHK H I T , CONTA M, EF IE LD , MA I N, RNOP OS , SH AOE * VNORML. 

COMMON /BLCK1S/ isectp 

C COMMON 3L0K1S IN D AT A IN , EF IEuD , E3U RF , MA I N, RND POS , S CPLOT * SHADE* ANO CONTAI 

COMMON /ELCK1S/ SP 1 ,SP 2 * SP 7 , SF 4 

C COMMON 8LCK21 IN CHK HIT , CO NTA M, EF IEL 0 . MA I N, RNDPOS . SHAOE » SPHIT. SCPL 

COMMON /ELCK21/ Z 1 » 2 2 » Z 2 « Z 4 » Z 5 

C COMMON BL0K2S CONNECTS CONTAM. EFIELD, MAIN AND TRAJEC. 
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COMMON /BL GK25 / EVCC<3» 

CCMHON BLOK 31 CONNECTS CONTAM, E'lELD AND MAIN. 

COMMON /BLCK31/ NRMISS 

EQUIVALENCE CISECTR.ISE 1 , < I P ANEL . IP J • I ILOST# IL I 

EQUIVALENCE (NEWPOSCl J .PCS NEW <1 H , (NEWEL tl »*VELNEV (1 ) > 

EQUIVALENCE tv EL VEC ( 1) . VELOLCt 1) ) ,< POSVECI1 I fPOSOLDt 1 ) ) 

EQUIVALENCE (VELMAG.VMAGl . (GMASS.GKAS ) , ( IS H ACE t IS H ) 

KCOUNTrC 

FATHLC r C. 2 INITIALIZE CONTAM PATH LENGTH ACCUMULATOR. 

IL^O 2 INITIALIZE ILOST REGISTER c OR NEW TRIAL. 

IHTtlCC 
THITt 2>=0 
NRMISSrO 
CELT AT;P. 

S Z rPOSOLD ( 3) 

CHGOLD=PRTCHG 

CO TO <20C. 220,220,10. ISE 

2SE IF( VELOLOl 3 ) • GT GO TO 1939 3 PARTICLE LOST 

C-0 TO 22C 

IP IF( VELCLDt 3).LE.O.OI GO TO 1399 

GO TO 22C 
20P CONTINUE 

I F t Z » CT .2 1-»Z 5 J GO TO 250 3 1 DEBYE LENGTH ABOVE 

220 CALL TRAjEC 3 INCREMENT PARTICLE POSITION. 

ishcld^ishade 

I POLO -IP 
ISEOLD=ISE 
KOLO=KCOUNT 
FTHOLD^PATHLG 

305 CALL SHADE 3 SEE IF PARTICLE IS SHAOEO AT THE NEW POSITION. 

SUM zC • 

DC 10C1 KK=1»3 

1CC1 SUM:SUM*tP0S0LD(KK>-F0SNEWtKKlt*»2 

PATHLG=PATHLG»SQRT<SUM) 3 THE ACCUMULATED PATH LENGTH. 
KC0UNT=KC0UN’*1. 3 INCPEHENT STEP COUNTER. 

TOPATHLG.GE.PATHMX) GO TO 1993 
iFtKCOUNT .EQ.5CCI GC TC 199S 

00 230 J -1 « 2 3 DETERMINE I c COORDINATES NEEO TRANSL. 

23C iFtPCSNEWC J) .LT.C.Ct GO TC 225 
GO TO 226 

225 CALL TRANSL 2 TFANSLATE PARTICLE BACK TO FIRST XY QUADRANT. 

236 CONTINUE 

? 4 P TF< -SE .E9 .4) GC TC 100C 

T”t IP.EO.l.AND. ISE.LE .2.4ND« S 0SNEWI3) .GT .73) GO TO 241 

SUOROUTINE CHKHIT DETERMINES WHE T HER THE NEW POSITION HAS CONTACTED 
A SPACECRAF T CONE SECTION. 

CALL CHKHIT 


C WHEN THE FASTICLE TRAVEP3ES THE CONE REGION AT THE TIP LOCATION 

C < Z 2. SP4 ) - THE C ALSE VALUE CF A RECONTAM HIT IS REVERSED TO A SAFE 

C HIT. 

C 

IFf ISE .NE .3 .CR .ISE0LC.NE.2.0F.IHIT (2>.NE.l J GO TO 242 
IHITt 1 J=1 
I HIT (2>=0 
30 TO 2001 
C 

C THE PARTICLE IS IN THE REGION OF A SOLAR PANEL AND IS MONITOREO BY 

C SUBROUTINE SPHIT FOR SURFACE CONTACT - SAFE FIT IF OCCURS. 

o 

241 CALL SPHIT — . 

#8S:SS* 
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C WHEN THE P ARTICLE LEAVES THE TOF OF THE SOLAR PANEL ANO BECAUSE OF 

C THE CHARGE-ELECTRIC FIELD CONDITION IS SUCH THAT THE PARTICLE IS, 

C SUCKED EACH THROUGH THE FANEL AnC MOVING DOWNWARO A SAFE HIT IS 

C REGISTERED. 

C 

IFIKCOUNT.EG.I.AND.I? .S3 . 1 . A V D . IP OLD . EQ . 1 . ANO . I ( NEVPOSI 3I.LT. 
1ZA.ANC .VELVEC13) »GT , C , C ) .OR. (NEW PCS ( 3 J . GT .24 .AND . 

2VELVECI 3) .LT.C.O) >) IHIT(1>=1 

242 IF< IHIT(1> .EG.C.AND.IHIT(2) ,EO.C> GC TO 1000 3 TO NEXT ITERATION. 

C ON TO STATISTICS EVALUATION 

GC TO 2CC 1 3 STATISTICAL EVALUATION 

C 

C INCREMENT THE VELOCITY AND POSITION VECTOR COMPONENTS. 

C 

:rcr continue 

lClC DO 1C11 LJ=1*3 

VELOLC(LJIrNEWVELtLJ) 

1011 POSOLD(LJ>=POSNEWCLJ) 

GC TO 5 

1030 IL=1 3 THE EJECTA HAS LOST ON THIS T RIAL. 

C CNE MONTE CARLO PASS COMPLETEC. 

2001 RETURN 

END 3 CCNTAM 
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SUBROUTINE C0NTM2 


'.UORCUTINC C0NT**2 

c • . • . >:»:t ..ri 


•••♦ 3/12/74 •••• D. ED GAR Si 8I0NETICS •••• 


CCNTM2 IS THE VERSION OF CONTAM THAT PRINTS THE TRAJECTORY INFOR- 
MATION FROM THE INITIAL GRAIN REMOVAL FROM T ME RANDOMLY SELECTED 
*0 CEL SURFACE -OCATICN TO THE DECONTAMINATION CONTACT. THE CHARGE 
AND THE ELECTPIC FIELD WRITTEN ON EACH LINE TELEOLOGICALLY CORRESPCNC 
TO THE PARTICLES EXPERIENCE AT THE POSITION COORDINATES PRINTEO. 

TK PATH LENCTH IS ACCUMULATED FOP ESCAPE CRITERIA EVALUATION ANO 
THE TIME OF FLIGHT INCREMENT FROM THE PREVIOUS POSITION IS INDICATED. 
T>-Is VALUE AFPPECI APLY CHANGES WHEN THE GRAIN TRAVERSES A SHAOE/SUN 
3CUN0RY. AN INDICATION OF THE T HPEE r OLD BALANCING MECHANISM INHERANT 
WITHIN ’HE SUBROUTINE TRAJEC CALLED HEREIN. THE VISUALLY INFORM- 
ATIVE COLUMN SH-P-SE REPRESENTS THE SHADE i PANEL AND SECTOR 
INDICTEE YIELDING QUICK IDENTIFICATION OF THE PARTICLE ENVIRONMENT. 

THE SUBROUTINE COORDINATES THE TRAJECTORY ANALYSIS FOR THE EJECTED 
CRAIN UNTIL THE GRAIN ACHIEVES ONE OF THE FOLLOWING CONDITIONS? 

1- THE GCAIN RECONTACTS AN ACCEPTABLE PORTION OF THE SPACECRAFT 

whereupon the variable tutiij is set to i - safe hit. 

2- THE GRAIN TRAJECTORY SURFASEES THE PRESCRIBED MAXIMUM PATHLENGTM 
OF PATHMX (20 «ETERS) - ESCAPEC. 


3- THE NUMBER OF TRAJECTORY INCREMENTS EXCEEDS 500 STEPS - ESCAPE. 

A- THE O'* A IN HAS RECONTACTED AN AREA OF THE SPACECRAFT CONSIDERED 
STEC ILE WHEREUPON IHITO is SET TO 1 - RECONTAM HIT. 


CURING THE EXECUTION THE PAPTICLE IS KINEMATICALLY PROJECTED FROM 
AN INITIAL POSITION THRU INCH EMENTEO POSITIONS USING THE PARTICLE 
ChARGEt THE ELECTRIC AND INERTIAL FORCES EXPERIENCED BY THE 
PARTICLE AT THE INITIAL LOCATION OF EACH INCREMENT. THE NEW POSI- 
TION IS SYMMETRICALLY MAINTAINED IN THE FIRST OCTANT AND MONITOREO 
FOR THE HIT CONDITIONS MENTIONED. 


C 

C 


THE INPUT PARAMETERS ARE 


P CS DL D • 3 1 
VELCLOl 3) 
CMAS 
VM AG 

ICE 

IP 


:sh 


POSITION VECTOR (INITIAL CONDITIONS). 
VELOCITY VECTOR (INITIAL CONDITIONS). 
The PAPTICLE (EJECTED MASS) 

THE VELOCITY MAGNITUDE 
DEFINES THE SE CT QP CF CONCERN. 

defines condition where: 


INPUT 

INPUT 

INPUT 

INPUT 

INPUT 

INTERNAL/OUT 


1 = PARTICLE UNDER A SOLAR PANEL. 

0 = PARTICLE NOT UNDER A SOLAR PANEL. 

CEFINEE CONDITION WHERE? INTERNAL/CUT 

1 ; PARTICLE SHADED FROM SUNtOUTSIDE WAKE: 
C = PARTICLE NOT SHADED FROM SUN. 

-1 - PARTICLE SHADED FROM SUN (INSIDE WAKE 

.«»« »*** 


REAL NEWPOSiNEWVEL 

DIMENSION PCECLD(I). PCSf.'CW(3). VELCLD(T>» VELNEW(3) 

COMMON GLCKOa CONNECTS CONTAM. MAIM. PARROT. SPHIT AND TRAJEC. 
COMMCN /?LCKr4/ CCEF ! 'F .CONST A. FRTCKC.FHIC 


C COMMON BLOKCB CONN EC’S CONTAM, DAT AI N » MAIN AND PARPOT. 

COMMON /BL0KC8/ DIM ,F AT H MX , S R , V E LM A C 

C COMMON 9L0KG3 IN CHK HIT . CO NT A M, MA IN, S HADE »S PH IT . T RA JE C ANO TRANSL. 

C 0 M MCN /SL0KP9/ NEW°CS(3) 

C COMMON 3L0KID CONNECTS CONTAM . MAIN, S PHTT .TRAJEC AND TRANSL. 

COMMON /eLOKIC/ NEWVEL(3). DELTAT 

C COMMON BLOK 11 IN CON T Am . EFI£;_ D . MAIN, P NOP OS ,S PHIT ,T RA JEC • V NORHL. 

COMMCN /BLOK 11 / P0SVEC(3) 
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C COMMON BLOK12 CONNECTS CONTAM. MAIN ANO TRAJEC. 

COMMON /3LOK1? / CMASS * VELVE C ( 3! .ILCST 
C COMMON BLOK 13 CONNECTS CHKHIT * CONTAM. MAIN AND SPHIT. 

COMMON /9LCK13/ IHITI21 

C COMMON BLOK11 CONNECTS CONTAM. EFIELD. MAIN. SHADE ANO TRAJEC. 

COMMON /BLOK1 1 */ IP ANEL . ISHACE 

C COMMON 3L0K15 IN CHK HIT • CONTA M. EF IELC . MA IN. RNDPOS .SHADE .VNORNL. 

CCMMCN /BLCK15/ ISECTP 

C COMMON 8LCK18 IN OAT AIN.EFIELD.ESURr.MATN. RNDPOS .SCPLOT. SHADE* ANO CONTAf 

CCMMCN /BLOK18/ SP 1 .SP 2 »SP 3 . SP «• 

C COMMON 3LCK21 IN CHK HIT . CO NTA M, EF IELC . MAIN. RNDPOS .SHADE. SPHIT. SCPL 

COMMON /BLOK 21 / Z 1 » 22 . Z 3 .Z«* . 15 

C COMMON BL0K25 CONNECTS CONTAM. EFIELD. MAIN AND TRAJEC. 

COMMON /BLCK25/ EVEC(3) 

C COMMON BLOK 31 CONNECTS CONTAM. EFIELD AND MAIN. 

COMMON /BLOK 31 / NRMISS 

EQUIVALENCE IISECTR.ISE ) . < I PANEL . IP ) . C ILOST. IL I 

EQUIVALENCE l NEWPOS ( 1 1 .PCS NEW I 1 > > » < NE V VEL 1 1 J » VELNEW U ) ) 

EQUIVALENCE < V C L VEC < 1 ) .VELOLDIll ) .(POSVECI1I .POSOLD(l) I 
EQUIVALENCE (VELMAG.VMAGI . ( CM AES. GMAS ) . (ISHADE.ISHJ 
1500 FCRMATC 1 X * IS* 2X » 5 I IP 01 2 • 5 * 2X) .12. IX. II . IX . II . 2X * 3« 1PG1D .*» 2X ). 

1 1 r C S . 3 1 

ISC1 F 0 RMATi 1 HP. * 4X * • ’’HE 3 ARTICLE ESCAPED’! 

1 S C2 FOPMATJIHI.?/. ’STEP PATH LENGTH X-POSITION Y -POSITION Z 
1-POSITION PART. CHG. SH-’ -S E ’ ,4X . ’E t X ) * * 8X * • ECY » • . 8X * 

2*E(Z)’»6X. *DT’ I 

1503 FORMAT! 1HC. ’FINAL STEP ’ . 2 X . ’ F INAL PATH LENGT H’ * 15X FINAL <X»Y*Z>* 

1.13X.* FINAL CHARGE’ »2X . * r I NA L-SH ACE- PANEL-SECT OR* * /» IX *16 *7X » 
21PG12.6.8X.M ’. 1PG10.4,’. ’.I’GIO.4.’, ’.1PG10.4. * > * * 2X » 1 PG12 .5* 

31CX.I2.5X.II.5X.I1> 

KC0UNT=0 3 INl T IALIZATTON OF THE WORKINC PARAMETERS. 

TL-C 

PATHLC " 

IHIT«1)=C 

TNIT(2!-0 

NR^ISSrC 

DELTAT=C. 

TIME :C • 

5 Z-P0S0LDI2) 

C HGOLD =FPT CHG 

GO TO C 20C.22C. 220.1DJ.ISE 

25C Ir»VELCLD«3l .CT.C.) CG TO 1959 3 PARTICLE LOST 

GO TC 22C 

1C Ir« VELCLO« 3) .LC.C.C) GC To 1995 

GO TO 22C 
2CC CONTINUE 

I FI Z. GT .21* ZS ) GO *0 2 5P 3 1 DEBYE LENGTH ABOVE 

:rc CALL TRAJEC £ INCREMENT PARTICLE POSITION. 

TT ME=TIVC*D"L TA~ 

TSHOLCrlSHAOE 
I P PLO-IP 
ISECLC=ISE 
KOLO=KCOUNT 
FTHOLO-FAOHLC 

3C5 CALL SHADE 0 SEE I c PARTICLE IS SHADED AT THE HEM POSITION. 

SUM-C. 

DO 1CG1 KK.= 1,3 

ICC1 SUM=SUM*IPOSCLO(KK)-POSNEW (KK) 1**2 

PATHLG=PA t HLG*S')RT<SUM) 3 THE ACCUMULATED PATH LENGTH. 

KMOD-MOCIKCCUNT.47) 

I F ( KMOD . EG. 0 1 HR ITE < o . 1 502 ! 

KC0UNT:KC0UNT*1. 3 INCREMENT STEP COUNTER. 

TFI PATHLG.GE.PATHMX) GO TO 1399 
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SUBROUTINE CONTM 2 


IF{KCeUNT.E3.5CCJ GO TC 1995 

00 230 J=1.2 a DETERMINE IF COORDINATES NEED TRANSL. 

27C IFJPOSNCWI JJ.LT.C.CJ GO TO 235 

GO TO 236 

235 CALL TRANSL 0 TRANSLATE PARTICLE BACK TO FIRST XV QUADRANT. 

235 WR ITEC 5. 1500) < 01 O.PT HO LOtPOS VEC t CH3 OLDt ISHOLD. IP0LD» ISEOLDtEVEC* 
1CELTAT 

240 IFtISE.EO.4) GO TO 1000 

1 Ft IP. EQ .1 .AKC .ISE .LE.2 • AND . F CSN'EU ( 3 ) .GT.Z3 ) GC TO 241 

SUBROUTINE CHKHIT DETERMINES WHETHER THE NEW POSITION MAS CONTACTED 
A SPACECRAFT CONE SECTION. 

CALL CHKHIT 

WHEN THE PARTICLE TRAVERSES THE CONE REGION AT THE TIP LOCATION 
(I2.sr4) - The FALSE VALUE of a recon t am hit is reversed TO a SAFE 

HIT. 

IF( ISE.NE. 3 .0 R. ISEOL C .NE . 2 1 GO TO 242 
D CLC=SQ°T tPOSOLOtlJ’RCSCLDtl )♦ POSOLC I ? » *POSOLD 12 )) 
RNEW=S3RTtP0SNEWtl)*P0SNEWIlJ *P OS NEW C 2) . POS NEW < 2 ) ) 

IFJRNEW.LT.SF4. AND. RCLD.LT.SF4) GO TC 2410 
GO TO 242 
41C I F IT ( 1 ) :0 
IHITJ 2 ) - 1 
GO TC 24 2 

THE FARTICLE IS IN THE REGION OF A SCLAR PANEL AND IS MONITOREO BY 
SUBROUTINE SP HIT FOR SURFACE CONTACT - SAFE HIT IF OCCURS. 

241 CALL SPHIT 

WHEN THE PARTICLE LEAVES THE T CP OF THE SOLAR PANEL AND BECAUSE OF 
THE CHARGE-ELECTRIC FIELD CONDITIONS IS SUCH THAT THE PARTICLE IS 
SUCKED BACK THROUGH 'HE PANEL ANO MOVING DOWNWARD A SAFE HIT IS 
REGISTERED . 

IFJKCCUNT .CO .1 .AND. IF .EG. I .AND .1 POL C .C Q.l . AND . It NEW P0SI3 » .L T. 

124. ANC.VFLVEC (3 ) .GT.2.C) .OR.t NE WPOS ( 3 ) . G T . 24. AND. 

bvelve: j *) .it .r .n > ) ihitiiui 

242 ITJ IHlTt 1) .23 .C.AND.IHITJ2) .EQ.O: GC TO 10C0 3 TO NEXT ITERATION. 

ON TC STATISTICS EVALUATION 

GO TO 2 r 0 1 3 STATISTTCAL EVALUATION 

INCREMENT THE VELOCIT* AND P CSITICN VEC'OR COMPONENTS. 

ICOr CONTINUE 
1C1C DC 1C11 L J -1 * 7 

VELOLOJ L J »=NEWVSL (LJ> 
iril PCSGLCILJ)=PCSNEWJL j) 

DC TO 3 

1599 IL -1 3 T HE EJCC T A WAS LC ST CN THIS TRIAL. 

W’lTEt St 1SCD) KOCO.PTHOLD.POS V2C. ' H3 OLD* ISHOLD. I POLO. ISEOLD.EVECt 
IDCLTAT 
WRITEfC.lSClJ 

WR ITE (E.15C3 » K COUNT .F AT HLG . AE W PCS . "T C H G » IS F ADE . IP ANEL * I SE C TR 
C ONE MONTE CARLO PASS COMPLETED. 

GO TO 2EC “ 

2CB1 CONTINUE 

IFJ IHITJ2) .EO.I) WRIT E ( 6 » 2CC4 ) 

MR IT El G» 152 3) K COUN T . P AT HLG . N EW PCS . P RT CHG . ISH ADE. IPANEL. ISECTR 
WOTTEJ6.15C4) TIME 


tfSKSSS 
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15C4 FC RM ATI 1 HC » * T HE PARTICLE TRAJECTOPY ELAPSEO TIME W AS * « G12. 3» 2Xt 
1 ’SECONDS * I 

FTA FORMAT! 1HC. *R ECONTAM IN ATICN HIT ON THIS TRIAL*! 

TOO RETURN 
END 
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sdbrout ine efielc a J. EAPENG0LT2. JPL »»lC/7 /74»» 

• EFIELD • 

CALCULATES THE VECTOR ELECTRIC? r IEL? AT A POSITION NEAR THE S/C 
USING THE SURFACE CONDITIONS ESTABLISHED 3Y ESL'RF. IN GENERAL. 

THE POSITION 13 TESTED FOR BEING NEAR A SURFACE IN WHICH CASE THE 
USUAL flat PLATE SOLUTION IS EMPLOYED WITH AN APPROPRIATE CEBYE 
LENGTH ANC SURFACE VALUE PER LENGTH AND SURFACE VALUE PER ESURF. 

THESE FACTORS DEPEND ON THE SUN/SHACC INCEX I ISHAOE ) ALSO. TWO 
EXCEPTIONS EXIS T t tlJ FOR POSITIONS NEAR ICOMPAREO TO APPROPRIATE 
CEBYE LENGTH) THE SOLAR PANEL, THE ALGEBRAIC SUM OF T HE SURFACE 
POTENTIAL OF THE ILLUMINATED AND DARK SIDES IS FORMED AND THEN 
THE FIELD IS CALCULATED ANALCGOUSLY. CD) FOR POSITIONS IN THE 
MAKE OF THE SPACECRAFT THE NEAR CASE AND FAR CASE ARE APPROXIMATED 
FOLLOWING AN ANALYSIS BY AL’FERT , ET AL.» SPACE SCIENCE REVIEWS. 

VOL. 2 (1353.) 

THE FAR FIELD CASE IS CALCULATED BY AN EQUIVALENT SPHERE CF AREA- 
AVERAGED POTENTIAL ANO THE DEBYE LENGTH 5ELCW THE SPACECRAFT OR 
AN AVER AGEC EFFECTIVE DEBYE LENGTH (SEE ESURF) ABOVE OR ALONGSIDE 
THE SPACECRAFT. 

OOU3LE PRECISION W. 30C2).3U2> 

DIMENSION COSJIA). SIN JC A I » TANJIA) 

COMMON BLOKU IN CON T AM • EFIELD. MAIN. RNOPOS »S PHIT »T RAJEC .VNORML. 

COMMON /BLOK11 / POSVECC3) 

COMMON BLCK 14 CONNECTS COSTAM. EFIELD. MAIN, SHADE AND TRAJEC. 

COMMON /8L0K14/ IP ANEL » ISHAC E 

COMMON BLOK1S IN CHKHIT.CONTAM.EFIELC.MAIN.RNOPOS.SHADE.VNORML. 

COMMON /BL0K15/ ISECTR 

COMMON 3LCK13 IN OAT AIN . EFIELD. ESURF . MAIN, RNOPOS »S CPLOT, SHADE. 

COMMON /BL0K13/ SP1 ,SP 2 »SP3, SP *1 

COMMON BLOK 13 CONNECTS DATAIN, EFIELD. ESURF. MAIN. RNOPOS AND SCPLOT. 
COMMON /BL0K19/ SPOl .SPC2.SP C3 

COMMON BLOK 20 IN CWHIT, EFIELD, MAIN, RNOPOS, SHADE AND VNORML. 
COMMON /SL0K20/ TANTH1 .TANTH2.T ANTH3 .T ANTH4 

COMMON BLOK 21 IN CHKHIT.CONTAM.EFIELO.MAIN. RNOPOS .SHADE. SPHIT, SCPL 
COMMON /8L0K21/ Z1.22.Z2.24.Z5 

COMMON BLOK 23 CONNECTS EFIELO. MAIN ANO SHADE. 

COMMON /BL0K23/ XYWAKE.2WAKE » JPCT 

COMMON BLOK 24 CONNECTS EFIELD. ESUR C- MAIN AND DATAIN. 

COMMON /BL0K24/AL AMBI2) .ANPHCT (2),AKPH0T(2 ). EFELC2 I ,EFED 12 ) .DEB. 

1AM ATI 9) »EL<2) .PHIEL 12 ) .PHIEOI 2) . ALAHAV. PHIAVE 
COMMON BL0K2S CONNECTS CONTAH* EFIELD. MAIN AND TRAJEC. 

COMMON /BLOK 2 5/ EVEC(3) 

COMMON BL0K2S CONNECTS DAT AIK. EFIELC .ESURF .MAIN AND PARPOT. 

COMMON /BL0K29/ AKTE. AKTP . ANE 

COMMON BL0K31 CONNECTS C? NT AM. EFIELD AND MAIN. 

COMMON /8L0K 31/ NRMISS 

EQUIVALENCE IX .POSVEC 111 I . I Y ,P OSVEC 12)) » IZ .POSVEC I 3 ) ) 

EQUIVALENCE (EX.EVECtll). IEY . EVEC 12) ) » (EZ.EVECI3) ) 

EQUIVALENCE IKSECTR. ISECTR) 

EQUIVALENCE IJ.JPOT) 

T AN J (1 ) =T ANTH3 
TANJI2)=TANTH2 
T AN J I 3 ) =T ANTH1 
TANJI4)=TANTH4 
CO 5C 1=1.4 

CCSJI I) =1./S3RT (l.«-TANJIII*TANJ(I)) 

5f SlNJIDrABSlTANJCI) )*COSJII) 

R= SQRTC X *X. Y* Y) 

J = 1 

EXrC.C 
E Y = 0 .0 

cz=c.c 
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PI-3.14159 

T F C IS EOT P-2) 100* 23G • 33C 
C SECTOR T IN CLOSE 

ICC IFJIPANEL.NE.il 30 T 0 110 

IFCAMATC3) .E3.ELC2) I J=2 3 AECVE SOLAR PANEL. 

DTST:A3SCZ-74 I 

IFCDIST »GE .ALAM3CJI) GO TO 1 1C 

E2 =C EFIL C J )♦ EFED Cl I * CE3/ ALA M3 IJ I) * EX P t- C IST/ A LA M3 C J I) 

RETURN 

11C IFIZ.LE.ZSI GO TO 120 
IFCR.GE .S D CI 1 GO TO 120 
C A30VE TOP CIRCLE 

iFCAMATtl) .EC.ELC2I I J=2 
DTST:A3SCZ-Z5) 

IFCOIST.GE.ALAMSCJI I GO TO 12C 
£Z=EFELC J)«EXPC-DIST/ALAM3C J) ) 

RETURN 

12C IFCZ.GE.Z51 GO TO 400 

C ALONGSIDE CONICAL SECTION 

IF C AMATC 2I.EQ.ELC2)) J-2 

OIST=C CR-SP2I* AES I T AN Jill ) *2- Z4 I * C OS J Cl) * 

IFCOIST.GE.ALAMBCJI) GO TO 430 

EMAG=EFELC J)*EXP C-OIST / ALAMB CJI I 

EX-EMAG»SINJC1)*X/R 

EYrEMA3*SINJ Cl >*Y/R 

EZ-EMAG»C0SJC1I 

RETURN 

20C CONTINUE 

IFCR.GT.SP3) GO TO 4C0 
IFCIPANEL.NE.1) GO TO 210 
DIST=ABSCZ-Z4I 

IFCDIST.GE.CZ4-Z3II GO TO 210 
IFC AMATC3I .EG.ELC2) I J=2 
C SECTOR II IN CLOSE 

EZ=-IEFEDC1I -tEFELC J» * A LA KB C J l/DEB I * EXP C-DIST/DEB I 8 BELOU SOLAR PANEL. 
RETURN 

21C IFCZ.LE.Z3I GO TG 220 

C NOT UNDER SOLAR PANEL BUT Z .3T. Z3 . AND .LE .SP3» CRUDE FIX WITH E=0. 

RETURN 

220 CONTINUE 8 NOT NECESSARILY UNDER PANEL. Z.LE.Z3 

IFCAHATC7) .EG.ELC2I I J=2 

OIST=CC R-SP4) • ABSCTAN JC 2)1 ♦Z-Z2I * COS JC2) 

IFCDIST ,GE .ALAMB C J 1 1 GO TO 4 CD 
C ALONGSIDE CONICAL SECTION 

ENAG=EFELC J l*EXP C-DIST/ A LA KB CJI) 

EX=EMAG*SINJC 21 »X/R 
EY=EMAS*SINJC2I»Y/R 
EZ-EMAG*C0SJC2I 
RETURN 

3CC CONTINUE 

C SECTOR III IN CLOSE TO WAKE EUT NOT INSIDE. 

ircZ.LE.O.) GO TO 400 8 BELOW WAKE. 

I F» ISHADE .E3 .-II GO TO 320 
RE = S3RTC X*X«-Y*Y»CZ-Z4)«CZ-Z4) I 
IFCRE.GE.ALAMAVI GO to 4C1 
DISTrCC P-XYWAKE)*ABSCTAN JC4I) ♦24-ZI* C0SJC4) 

IFCDIST. GE.ALAMAVI GO TO 4CC 

EMAG-EFEDC11»EXPC-0IST/DE3) 3 LEADING COEF. SHOULD BE HAKE SURF. 

EX-EKAG*SINJC4 l«X/R 
EY -ENAG*SINJC 41 »Y/R 
EZ=-EMAG*C0SJC4> 

RETURN 

32C CONTINUE 3 INSIDE WAKE. 

TFCZ.LE.Z1) GO TO 330 
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C NEAR BOTTOM CONE. CRITICAL AREA 

DIST=C < P-SP4)» ABS (TAN J( 3 ) ) ♦! 2- Z ) • CDS J( 3) 
EMAGzEFEDCIMEXPI-DIST/DEB) 

EX=£MAG*SINJC 3) *X/R 
EY=EMAG*SINJI3)*Y/R 
EZ=-EMAG*COSJ 13) 

IFIDIST.LT.O.C1) NRMISS =1 
RETURN 

33C CONTINUE a BELOW S/C IN WAKE. 

FIELD PER J3B 9/10/74 FOR CYLINDER. 

VALID CLOSE TO S/C BOTTOM ONLY. 
ACYL=(Z1-Z4»ZWAKE)*XYWAKE/ZWAKE 
DISTr ABS CZ1-ZI 
W-2.4*R/ACYL 

IFCCIST.GE ,ACYL» GO TO 39C 
CALL BJYOK W* 1. 1* 80. Bl) 

EMAGR0=2.4»EFED(1>* (DEB/ ACYl )• EX PC-U ) 
EZ=-EMAGR0*B0<1) 

IFCR.LE.O.I RETURN 
EX=EMAGR0»X/R*B1«1) 

EY=EHAGR0*Y/R*B1(1I 
NRMISS=1 
RETURN 

390 CONTINUE a In WAKE but far from bottom. 

ehagro=-akte/deb 

EX=EMAGRO*X/R 
EY=EHAGRO*Y/R 

EZ=I AKTE/ACYL 1* ALOG(OEB/ ACYL) *EXP(-W) 

RETURN 

400 CONTINUE 

C EQUIVALENT SPHERE 

RE^SQRTI X*X* Y *Y*I Z-Z4I * ( Z-Z4) ) 

401 ALAM=AL AHAV 

EMAG=PHIAVE*(1./RE)*(1./RE»1./ALAH ) *XYW AKE 
IFIRE. LE.XYWAKE) GO TO 470 
EMAG=EMAG*EXP I- (RE-XYWAKEJ/AL AM) 

470 CONTINUE 

EX=EMAG*X/RE 

EY=EMAG*Y/R 

EZ=ENAGMZ-Z4I/RE 

RETURN 

END 
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SUBROUTINE ESURF 


SUBROUTINE ESURF 3 J. 8ARENG0LTZ, JPL»*10/7/74 ** 


» ESURF • 

CALCULATES SUNLIT SURFACE POTENTIAL PHIEL FOR AL AND SI02 AND THE 
CORRESPONDING SURFACE FIELD EFEL AND THE DISTANCE PARAMETER ALAMB 
PER JPL GTP. VOL. 3 NC It APRIL 13?3t AND SETS SHADED SURFACE 
POTENTIAL tPHIEDl TO -3*ELECTR ON TEMPERATURE (IN VOLTSJ AND 
SURFACE FIELD tEFEDI TO THE FOTENTIAL CIVIOED BY THE DEBYE LENGTH 
(DE3.1 THIS ROUTINE ALSO CALCULATES . AREA- AVERAGED SURFACE POTENT- 
IAL (PHIAVE) AND AN EFFECTIVE DEBYE LENGTH t ALAMAV • I 

DIMENSION Af:l CN C 2 * »AA(S) 

COMMON 3L0KO3 CONNECTS ESURF, MAIN ANO RNDPOS. 

COMMON /BLOKC3/ ASP » ATP » ATS , PR It AT OT . ABP t ABCS • ATI 

COMMON BLOK 13 IN OAT A IN t EPIELDt ESURF , MAI Nt RNDPOS t S CPL OT t SHADE . 

COMMON /BLOK 13/ SP 1 t SP 2 t SP 3 t SP 4 

COMMON BLOK 13 CONNECTS DATAINt EFIELDt ESURF t MAINt RNDPOS AND SCPLOT. 
COMMON /BLOK 13/ SPC1 tSPC2tSP C3 

COMMON 3L0K24 CONNECTS EFIELD , ESURF MAIN AND DAT AIN. 

COMMON /BLOK 24 /ALAMB (2 J t ANPH CT ( 2 1 1 AK PHOT C 2 It EFEL 1 2 > t EF EO C 2 ) tOEB t 
1 AM ATI 9 J t EL(2» tPHIEL <2 ) ,PHIED(2» t ALAMAVt PHIAVE 
COMMON BL0K2S CONNECTS C AT AI Nt EF IELC t ESURF t M A3 N AND PARPOT. 

COMMON / BLOK 2 9/ AKTE t AKTP t ANE 
PI=3. 14159 
DO IOC I-lt 2 

ri=akte/akphot ti j 

R2=ANPH0Tt D/ANE 

IFCR1.LE.1.0.0P.R2.LE.R1) GO TO 900 
PHIEUIJ=AKPHOT II) • < ALOG <R2 )- 0.5» AL0GCR1) ) 

AN TON ( 11 = t ANE/2.1* tSGRT tRIJ ♦ 1.CM1.E-6 

CO=l ANPHOT tll/SRRTfAKPHOTtU) - ANE/ SORT ( A KT E 1 1 *l.E-6 

C1=-1.5*S0RT (PIJ*ANICN(I)/CP 

C00-8./3.»S3RT<PIltCC*1.44E-7 

C2 = SQRTfC0C/2. 1/2. 

Cl AND C2 AND FORMS OP PHI AND E TAKEN FROM BARENGOLTZ AND BAUERLEt 
JPL 3TR*LYt VOL. 3t NO 1 C1973> EXCEPT EX100 FOR V/M. 

C4=S3RTIPHIEL (I) J 
C5=SGRTf C4*C1J 
XMAX=C5/C2 
XDELT=XMAX/irC. 

EFELt I>=4C0.* CD *C4*C5 
XrXDELT 

DO 200 J-l t 1 00 

EEX=4CC .»C2» (-C1 •t(C5-C2«X J •• 2 . )• IC5-C2»XJ 
RT EST=EEX/EFELl I) 

IFCRTEST .LE.C.3G79) GO TO 3CC 
XrX*XOELT 
200 CONTINUE 

WRITE1 G t 902 ) 

STOP 

300 AL AMB( I) -X/100. 3 CONVERT TO METERS. 

PHIE0II»=-3.«AKTE 3 SHADE SURFACE. 

EFEOt I»=PHIEDIIJ/DEB 3 SHADE SURFACE. 

IOC CONTINUE 
AAt 1):ATP 
AA(2t=ATS 
AAt 3>:ASP/2. 

AAt 4 JrAAtD) 

AAt S):PI*SP2*SF2 
AAt G)=PI*SPC2*SPC2 
A A t 7 ) : AI I 
AAt 8):ABCS 
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AA(R>=A2P 

C FAR FIELD CADE • ALL SECTORS (APPROX.) 

ALAMAV=O.C 

PKiAVE=r.r 

ASC-C.C 

DO 44G 1=1*3 

J-l 

IFt AM AT (I) .EG. EL 12)1 J = 2 
PHIAVE=PHIAVE 4 A A ( I)*PHIEL ( J) 

ASC = AS C 4 A A (I ) 

44 r ALAMAV-ALAMAV 4 ALAM3( J)*AA(I) 

DO 450 1=4*6 

PHIAVE=PHIAVE4AA(I)*PHIED(1) 

ASC = AS C 4 A A (II 

45D ALAMAV = ALAM4V 4 DEB*AA(I) 

J-l 

IF( AM A T < 7)«EC«EL(2)9 J = 2 
PHIAVE=PHIAVE 4 AA (7I*FHIEL( Jl 
AL AMAV=ALAMAV 4 AA (71 *ALAMS ( J) 

A SC = AS C 4 A A (7) 

DO 4 SO 1=3*9 

PHIAVE=PHIAVE 4 AA (I)*PHISD(1) 

ASC=ASC 4 AA(I) 

46 T ALAMAV=ALAMAV4AA(II*DEB 
PHIAVF=PHIAVE/A$C 
ALAMAV=ALAMAV/ASC 

WR IT E( S * 10001 AL AM 9* EFEL* EFE3 * DE9* PMIEL* PHIEO * A LAMAV * PHIAyE* 

1AMAT .EL 

100C r 0 RM A T t 4X * * * • ALAMB • **tlOX*** • EF EL • * , *10X* 4 * » EFEO • 

12X. 4 DEB*./*1X* 7I1PG9. 3*2X1 *//.4X*** • FUEL • * f *9X*** • PHIEO • • 
2 • 

3*6X* 4 * AL AM AV • • PHIAVE * • * ✓ * IX * 8 ( 1 PG9 .3 * 2X I*// *4X * *♦ • AMAT ARRA 
4Y • * , *19X* t * • AMAT ARRAY • • * »/ . IX * 9( AC* 2X ) * / / * IX* • • • EL ARRAY 
5 • * *•/ *XX 1 2 ( AS* 2X I I 
RETURN 

9CC WRITECS.9C1I EL (I I 
STOP 

9C1 FORM AT (IX. * p HO TO ELECTRIC DAT/ INCONSISTENT WITH PLASMA DATA* SEE M 
IAIN FOR ’.AS) 

902 FORMAT ( IX * * EFFECTIVE SUNLIT CEBYE length; COULD NOT BE FOUND BY ESU 
1RF • ) 

ENO 
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MAIN DRIVER (VERSION - FULBPS) 


C •••* NASA JFL **•» D. EDGARS IBIONETICS). J. BARENCOLTZ (JPL) 

Q ***** ***** ***** 

REAL NU ♦ NEWPCr* N EWVC L» IVESC 


DIMENSION UH€J C5 * * WNE T 2 ( 5 * • V 2 ( 36 I * IS FUN (36 ) *W ANGDS I 1 8 } * 

1 V ANGDS t 13 ) 

DIMENSION AN 6DE6 (18) 

DIMENSION P0SHIDC3) 

DIMENSION V AP (35 ) * PKACC<3f> 

DIMENSION ANG DIS ( 1 8) ♦ FLUXMKItlO) 

C ***** ***** ***** 

C 

C this Is the DRIVER PPCGRA* fCR THE BYPASS MOCE OF THE SPACECRAFT 

C RECONTAMINATION MODEL. THIS MODEL HAS BEEN BRIEFLY DISCUSSED IN 

C THE JET PROPULSION LABOPATORY DOCUMENT 90C-675. SECTION III. 

C SEPTEMBER. 137*. 

C 

C THE DISTRIBUTIONS OF THE MASS AND VELOCITY FOR THE METEOROIOS IN 

C THE SPATIAL REGION BETWEEN THE EARTH ANC MARS HAVE BEEN C A TAGORIZED 

C INTO TEN MASS AND THREE VELOCITY GROUPS. .SE PE RATE EXECUTION OF 

C THE SURFACE EFFECT MODEL { Y A KG 1/THlNPL > HAS PRODUCED DATA ENSEMBLES 

C WHICH THIS PRESENT VERSION READS AS INPUT. THE PARTICLE RELEASE 

C MODEL tRELEAS) CONTAINED HEREIN COMPUTES THE STATISTICAL QUANTITIES 

C FOR PARTICLE RELEASE PROBABILITIES. 

C 

C THE METEOROID MASS AnD VELOCITY GROUP UNDER CONSIDERATION THEN HAS 

C RANDOM POSITIONS FOR IMPACT DETERMINED. AT EACH OF WHICH RELEASED 

C GRAINS ARE KINEMATICALLY FOLLOWED ON THEIR TRAJECTORIES. THESE 

C TRIAL GRAIN HISTORIES lEAD Tc ONE OF SEVERAL SITUATIONS* THE GRAIN 

C CAN RECONTACT THE SPACECRAFT ON A REGION WHERE THE STERILITY IS OF 

C NO CONCERN. THIS IS REFERRED To as a safe-hit. the grain may 
C ALSO ESCAPE THE VICINITY OF THE SPACECRAFT WHEREUPON THE NUMBER OF 

C AND DIRECTION of the escape IS NCTEC for later TABULATION. the 

C IMPORTANT OATA GENERATED WITHIN INDICATING THE NUMBER OF GRAINS 

C REACHING THE PRESCRIBED STERILE AREAS ARE Called RtCowT. HITS. 

C 

c ••».» ***** •*••• 

C 

C THE COMMON BLOCKS LISTED INTERCOMMUNICATE ESSENTIAL DATA AMONG THE 

C SUBROUTINES 

C 

C COMMON BLOK 01 CONNECTS DATAIN. MAIN AND RELEAS. 

COMMON /BLOKGI/ AKM »DRHO . K1 » FR BV EL f 35 . 1C ) » CRI1C I.SIG.DDIA ClC I 
l.ANPART* ANN0RM.PR80IAI10I 

C COMMON BLOKOZ CONNECTS CHKHIT. MAIN. RNOPCS AND SHADE. 

COMMON /BL0K02/ ALU) 

C COMMON BL0KC3 CONNECTS ESURF. MAIN AND RNOPOS. 

COMMON / BLOK 0 3/ ASP . ATP . ATS .P R I . ATOT . ABP. ABCS . All 
C COMMON BLOKO* CONNECTS CONTAK. MAIN. PARPOTt SPHIT AND TRAJEC. 

COMMON /BLOKO*/ CO EFRP .CONSTA • PRTCH3 . PHIO 
C COMMON eL0K05 CONNECTS OAT AIN. MAIN* RELEAS AND YANG1 

COMMON / BLOK 05/ KO.PM ASS.PRHO .RSI 35) *FVAPC 35) »FPKACC« 35 > • RHO 
C COMMON BLOKOG CONNECTS DATAIN. MAIN AND YAN61 . 

COMMON /BLOKOG/ E. H.NU.PO * ISKIP 
C COMMON BLOKCe CONNECTS CONTAK. DATAIN. MAIN ANO PARPOT. 

COMMON /BL0KC8/ DIM »P ATHMX. SR .VELMAS 

C COMMON ELOKC? IN CHKHIT .CO NT AM.MAIN.SHACE.SFHIT.TRAJEC ANC TRANSL. 

COMMON / BLOK C 3/ NEWP0SI3) 

C COMMON BLOK 1C CONNECTS CONT A M» M A IN . SPHIT »T R A JEC ANO TRANSL. 

COMMON / BLOK 10/ NEWVEL (3) . DILTAT 

C COMMON BL0K11 IN CONT AM . EF IE LD » M AIN , RND POS .S PHIT . T R A JEC » V NO RML . 

COMMON / 3L0K 1 1/ PCSVECI3) 

C COMMON BLCK12 CONNECTS CONTAK. MAIN AND TRAJEC. 

COMMON / 9L0K 1 2/ GMASS .VELVECI3) .ILOST 
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c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


c 

c 

c 

c 

c 


COMMON BLCK13 CONNECTS CHKHIT. CONTAM. MAIN AND SPHIT. 

COMMON /BL0K1 3/ IHIT(2) 

COMMON BL0K1 4 CONNECTS CONTAM. EFIELC. MAIN* SHADE AND TRAJEC. 

COMMON /3L0K14/ IPANFL * ISHADE 

COMMON BLOK15 IN CHKHIT , CO NT A M.E PI ELC . M AIN . RNDPCS ,S H ADE . VNORML . 

COMMON /BLOK 15/ ISECTR 

COMMON BLCK16 CONNECTS MAIN* RNDPCS »SHADE AND VNORML. 

COMMON /BL0K16/ NPNL 

COMMON BLOK17 CONNECTS CAT AIN* MAIN AND Y ANSI * 

COMMON /BL0K17/ PVEl*ACCMIN 

COMMON BLOK16 IN DAT AIN »EFIE LC *ESU RF .MAIN . RNDPCS » SC PLO T* SHADE . 

COMMON /BLOK1 8/ SP1.SP2.SP3.SP4 

COMMON BLOK19 CONNECTS DATAIN* EFIELC* ESURF.MAIN. RNDPOS AND SCPLCT. 
COMMON /BLOK 19/ SPOl . SP02 »SPQ 3 

COMMON BL0K2C IN CHKHIT* EFZELD* MAIN* RNDPOS* SHADE AND VNORML. 
COMMON /BL0K20/ T ANTH1 *TANTH2 »T ANTH3 * T ANTH4 

COMMON BL0K21 IN CHKHIT . CO NT AM.EFIELD* MAIN* RNDPOS *S HADE. SPHIT* SCPL 

COMMON /BL0K21/ Zl. 72 * 23* Z4 »Z 5 

COMMON BL0K22 CONNECTS MAIN AND VNORML. 

COMMON /3L0K22/ VELNRM 13) 

COMMON BL0K23 CONNECTS EFIELC* MAIN AND SHADE. 

COMMON /BLOK 2 3/ XYWAKE.ZWAKE » JPOT 

COMMON BL0K24 CONNECTS EFIELC* ESURF AND MAIN. 

COMMON /BLOK 2 4/ AL AMBC2) • ANPH0TC2I .AKPHOT (2) .EFEH21. EFEOC2U OEB* 

1 AMAT ( 9 ) »EL (2 ). PHIEL (2) *PHIED 12 I * AL AMAV * PHIAV E 
COMMON BLOK 25 CONNECTS CONTAM* EFIELD* MAIN AND TRAJEC. 

COMMON /BL0K25/ EVEC(3> 

COMMON 3L0K2S CONNECTS DATAIN* MAIN AND PARPOT. 

COMMON /BL0K26/ MAT * AVP 

COMMON BLOK 29 CONNECTS DATAIN .EFIELD* ESU RF, MAIN AND PARPOT. 

COMMON /BL0K29/ AKTE.AKTP* ANE 

COMMON BLOK 30 CONNECTS DATAIN AND MAIN. 

COMMON /BL0K3C/ADRH0(3>*D1»D2*D3*D4*N1*N2*SRMAX.SRMIN 
l.EPSILN . ELECT ♦ AP VEL < 3) .TMISS* KPLOT 
COMMON BLOK 31 CONNECTS CONTAM* EFIELD AND MAIN. 

COMMON /BL0K31/ NRMISS 

EQUIVALENCE IFVAP(l>»VAPCin .(FPKACC(l).PKACC(l l> 8 MAIN* RELEAS 
EQUIVALENCE (K1D.JSIZE) » (KOD* JVEL) 

PI=3. 14159 


the next approximately fifty lines of code calculate the necessary 

MODEL WORKING PARAMETERS FROM THE USER SPECIFIED SPACECRAFT 
INFORMATION ENTERED IN THE BLOCK DATA ELEMENT DATAIN. 

DRHO=ADRHO (MAT > 

ALPHAS AT AN (2.*1 »6E4*SGRT (AKTP •/ IAVP*S QRT (PI) I ) 8 WAKE HALF ANG. 

ORPIC6=DRHO*PI/G. 

WRITE(S.l) 

1 FORM AT ( 1H1 . *THE FOLLOWING ARE GENERAL SPACECRAFT MODEL DATA*) 

TANTH1=D1/(SP4-SPC3» 

T ANT H2 s-D 2/ (SP4-SPC2) 

TANTH3=-D4/ (SP2-SP01) 

T ANT Hfl :T AN (P 1/ 2 .- ALPHA I 

WRITE( S.7CD) T ANT HI • T ANT H2 * TA NT H3 *T ANTH4 
7CC FcRMAT(lHr,*THE TANGENTS FOR CONE ANGLES A RE * . 4 ( G9 . 4 . 2 X ) ) 

C 

XYWAKE = (SF4*SQRT «SP4*SP4*8.*SP1* (SP2-SP4I/PI > >/2 . 
ASP=1G.*SP1*(SP3-SP2) 

ATP=PI*<;fc1»«2 

ATS=PI*( SPC1*SP 2) *S3RT( 04* D4MSP2-S PCD * ISP2-SPC1J) 

AI = ASP-*ATS*ATP 

ATI=PI.( SP02*-SP4)*SGRT(D2»D2* (SP4-SP0 2) * (SP4-SP02) ) 

AT0T=AI*AII 
PRIsAII/ATOT 
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ABP=PI*SPC3**2 

A3CS=PI*C SP03-*SP4)*SQRTCD1*C1*CSP4-SP03) *CSP4-SP03)) 

VRITE<G.702) ASP. ATP.ATS.AI. AII.ATCT.PRI.ABP. ABCS 
7C2 FORM ATI 1 HC. • ASP = * * G7 . 3. 3X . ’ ATP = * . G7. 3 . 3X . » ATS = • . G7 .3. 3X . • A I=* . 

lG7.3.3X»*AIIr*»G7.3»3X. * ATOT r * ,G7 . 3 . 3X » *PRI= • *G7 .3 » 3X ♦ ’AB P= • » G7. 3 
2.3X.*ABCS=’.G7.3) 

SQTARG=EPSILN*AKTE/ CANE* ELECT) 

SUMFL = 0.0 
IVESC = O.C 
SMVESC=0. 

ACCUMlrC. 

ACCUM2=0. 

ACCUM3:C. S FOR NRMISS. SEE EFIElD 

DO 3 LM=1»18 

3 ANGDISCLM):0. S INITIALIZE THIS ARRAY EXPLICITLY. 

10 SR =10. * * < C. 5* CALCGICCSRMlN) *ALOG10CS RHAXI ) ) 

Z1=SQRTCSQTARG)*SR 3 DEBYE LENGTH AT RADIUS SR 

Z2=Z1*D1 

Z3=Z2*D2 

Z4=Z3*D3 

ZS=Z4*D4 

DEB=Z1 

HRlTECG.72CIZl.Z2.Z3.Z4.Z5 

720 FORMAT! 1 HO. *Z 1=’ . G9. 4. 2X . ’Z2= ’. G9 .4, 2X. »Z3=’. G9.4.2X. ’Z4 = ’» G9.4» 

12X » *Z5r *» G9.4I 

ZH AKE=Z4— XYVIAKE*T AMTH4 
URITEC6.7C1I XYUAKE.ZUAKE 
701 FORMAT f 1H0. • XYH AK E= • • F9 .4. 4X . *ZU AK£= F9.4 > 

C THE al*s are The z-axis intercept for the appropriate cone. 

ALC 1)=Z2-SP4*TANTH1 
ALC2)=Z2-SP4*TANTH2 
ALC 3)=Z4-SP2*TANTH3 
ALC4)=ZHAKE 
HRITEC6.721IAL 

721 FORM AT 1 1HC. ’THE AL ARRAY CO NT AINS ’ . 4 C69 .4, 2X ) ) 

CALL E3URF 

4 READ C 5.1111 ) N2.KFL0T 
1111 FORMATC 2 110 ) 

5 READ C 5.1600. ERR=50.END=5G) I HCU P» KO.IGML. IGMU . PVEL . PMA SS.IM V.T**“ 

1SM1 

SMTESTrO. 

HRITEC6. 1609) IMM.TMV 

1 609 FORM ATI1H1, ’THE INP'JT METEOR CIO CASE DATA CARD TOR MASS GROUF’. 

114. 2X. ‘AND VELOC. GROUP’. 14) 

WR ITEI 6.1E1C JIHICUP. KO . IG ML. IC MU . PV EL. PMASS . IMV . IMM . SMI 
ICOr FORMAT! 41 I3.2X) . 10X . 2 E12 . 5 . IX . 2 C 12. EX ) . E9 .4 ) 

IE 1C FORM ATC1HG.4 Cl 3, EX > . 1CX . 2E 12 .5 . IX , 2 1 12 . 2X ) . EE .4 ) 

IFtKPLGT.CO.l » CALL SCP LO T C I GC CN) 3 SPACECRAFT PLO T . 

IF! IGCCN.EG.-l > STOP 
2C FL UXM1C TMV.IMM)=.25*PVEL*SM1/Z. 

Fly OT rFLUXMl (IMV .EMM ) • ATOT*T KTSS*8.64C4 

WRITEI5.7C5J IMM.IMV .PM ASS VEL. FLU X Ml C IMV . IMM ) . FLT 0' 

705 F0RMATC1H1. ’PASS GRCUF = ’ . 1 3 . EX , • VE L OC. GROUP= * .13 . //.IX, • THE METEO 
IPO ID MASS IS* »1 p SI0.4,2X»’THE METEOROID VELCC. IS’*1PG10.4. / 

2/. IX. ’WITH A GROUP FLUX 0F»»1FG1C.4. 3X.’AND’.1X, 

3*THE TOTAL METCOROID IMPACTS OF THIS TYPE EXPECTED ARE’ . IX . G3.4 J 
WPITE (6.1599) 

READ! 5.1EC1)ID’JMY1.RS Cl) .PXACCCI) .VAPCI) 

WR IT EC6. 16CDIDUMY1.PS Cl) .PKACCC I) »VAF Cl) 

SIC CONTINUE 

15 99 FCRMAT(1HC»////»2X» ’I’.llX. »RC I) • . 9X • • PKACC C I ) • . 9X , »VAF C I ) ’ ) 

1601 F0RMATCI3.3X.1P3E15.5) 

TFCFLTCT.LT.r.Cl) GO TO 24C 

27 CALL RELEAS 3 DET . PA RT. RELEASE PARAM. 
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CO 8GG1 1=1*18 
3PC1 VANSDS! 11= O.C 
DO 8CCC 1=1.5 
3C0C WNET2!I»= C.C 
C ALGORITHM TO PACK VAT ARRAY 

VAP2!1>= ABS i VA° ! 1 ) ) 

X J=2 

TSFUN1 11=1 
DO 81C0 1=2. KC 

DO 8200 II=IJ.KO 

isfunii>=ii 

TF( ABS! VAP1II)) .LE. C .5* ABS (V AP ! I J-H > > GO TO 8300 
820C CONTINUE 
830C IFIX=ISFUNII) 

VAP2 ( I )= C.5*!A8S!VAP!IJ-1»»4ABSIVAP!IFIX) I) 

IQUIT=I 

IF1ISFUN1II ,E9. K0» GO TO 84 CC 
IJ=ISFUNlI)4l 
810C CONTINUE 
9400 CONTINUE 

IQUIT=I0UIT41 
VAP2I IQUIT»=A8S! VAP1K0T) 

EXIT KITH GRCUP POINTING ARRAY IN ISFUN1U .I8UIT IS NO. OF GROUPS* 

ANO VAP2tI) IS GROUP VEL. ARRAY 
DO 955 IT = 1.IHICUP 

955 OUNNY=R ANDNO ! 1. »0« ) 9 INCREMENT THE SYSTEM RANDOM NUMBER GEN. 

DO 4C J=1.N2 a INNER MONTE CARLO INDEX 

WRITE! 6.2) IMM.IMV.J 

2 FORM AT C 1H1 »/// »1X » ’MASS GROUF • » I3.5X • • V ELOC . GROUP* . 13 *5X . 

IMPOSITION NO.*. 15 I 
DO 4C00 1=1.5 

4000 WNET! 11=0.0 

DO 4001 1=1.18 

4001 WANGDSt 11=0.0 

CALL RNDPOS a DETERMINE IMPACT LOCATION 

CO 23 K=l,3 

rcsHic «ki zpo^vrc ck » 

29 NCWPOS! H )=P9E VEC('< ) 

CALL SHADE 
TSHHID=ISHA9E 
TPHIC = I r ANEL 
TSEHID= IS2CTR 

CALL VNORML 2 NORMAL UNIT VECTOR aPOSVEC 

WRITc(S.7n3) POSVEC. 'SHADE. IPA NIL. NPNL. IS ECTR 

7C9 F0R«AT(lHr, *THE PANDCM FCSITICN VECTOR IS*. 

13! G3. 4.2X) .//.IX. *ISHACE=**I2»2X» 

2*IPAnEL-*.I 2»2X. *NPNL=*«I2»?X» *ISECTRr*,l2 J 
WR ITE! S.71CI VELNRM 

71C FORMATCIHC. *’HE NORMAL VECTOR FOR THIS CASE IS * . 3 !G1C .4 »2X M 
IF! ISHADE.N'E. 01 GO TO 31 
CALL E FIELD 
PHID = PHIELI JPOT > 

GO TC 32 

31 PHI0=PHIED1U 

32 CONTINUE 
PHTDE=PHIO 

CO 5CC0 jsize=igml.igmu a jsize.eq.kid grain loop 
DIM=DCIA!K1D» *1 .GE-6 

C0EFRP=-7.E-G*PI».25*!DIM/SR>**2 a RADIATION PRESSURE FORCE IN SUN 

CONST A= 2 «*P I* EP SILN •DIM 3 CONVERT Ap«I TO AQ 

GMASS=DRPI0S*1DIM**3I 

DC 6000 JVEL=1.IQUIT 5JVEL ElUIV. KOD .GRAIN VEL LOOP 
C CALCULATE RELATIVE PROBABILITY OF JST2E.JVEL PAIR 

IFIJVEL .E8.1> 60 TO 9001 
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iFtJVEL ,EO. I3UIT) GO TO 9CC3 

IFIX=ISFUN(JV£L) 

iFIXlrlEFUN (JVEL-1) 

P3VElC:PR 80IA I J SIZE) • <PR B VEL( IF IX . JS I ZE) - PRBV EL t IFIX 1 » JSIZ E I) 

GC TC 9CC2 

9C01 PBVElC:PR3CIA(JSIZE)*PRBVEL(1. JSIZE) 

GO TC 9CC2 

9003 PBVELC: PRBDIAt JSIZE) • (l.-»RBVEL(KC» JSIZE) ) 

9002 CONTINUE 

WEIGHT:PBVELC*ANNORM 3T0TAL GRAINS THIS SIZE ANO VEL REMOVED 

V ELM AG: V AP2 ( JVEL ) 

DO 35 NN=1.3 3 VELOC. VECTOR (INITIALLY) 

POSVEC(NN):POSHID(NN) 

35 VEL VECI NN ) :VELM AG*VELNRM INN) 3 BUILD INITIAL VELOCITY VECTOR. 

PHIO:PHIDE 
PR TCHG: CONST A »P HIO 
ISHAOE:ISHHID 
IP ANEL: IPHID 
ISECTRiISEHID 
ATEST:WEIGHT*FLTOT 
SMTEST :SMTEST ♦ATEST 

IF( WEIGHT .LE. .001.0 R . ATEST. LE .0 .05 ) GO TO 5999 
IFC GMASS.LE.C. ) GO TO 800 3 DEBUGGING 

VMAX=SQRTI-1D.*C0EFRP*PATHMX/GMASS) 3 5 TIMES THE MINIMUM ESCAPE VEL. 

IFIVMAX.LT .VELMAG. AND. NPNL.NE.-l) GO TO 36 

iFCVELVECf 3>.LE.0..AND.NPNL.NE.-1> GO TO 310 3 DEBUGGING 

CALL CONTAM 3 TEST CASE TRAJECTORY ANALYSIS 

IFINRMISS.E3.1) WNE’ ( 3) = VNET ( 3 ) ♦ WEIGHT 

if( :lcst.ec.:j go tc 39 

TF(IHIT(1).E3.1) WNET ( 1 J :WNET til ♦WEIGHT 
IFC IHITtZI .E3.ll WNET (2) :WNET(.:>+U EIGHT 
IF( IHITt 2I.E3.U GO TO 30CC 
GO TC 6C0F 
3CCC CONTINUE 

C THIS PRINTS ’HE TRAJECTORY OF THE P RE DEEDING RECONTAMINATION EVENT. 

DO 3500 NN: 1 . 3 
FCSVECtNN):PCSHID«NN» 

35CC VELVECtNN):VELNAG*VELNRM(NN) 3 3UIL0 INITIAL VELOCITY VECTOR. 

P RIO :P HIDE 
P.TTCHG:CONST A *»HIO 
ISHADE:ISHHIC 
T p ANEL:IPHIO 
ISECTR:ISEHIC 
WRITE16.2CC0)J» JSIZE. JVEL 

2C0C FORM AT (1H1 » 'POSITION'.IE.AX. ’GRAIN SIZE GROUP* .I6.4X. 
l’GRAIN VELOC. GROUP*. 16) 

WRITE (6.7111 DIM.COEFRP.CONSTA.PRTCHG.GMASS 

l.PRBDIAt JSIZE ) 

711 FCRMaT( 1HC,///.1X. »DIM:».G9.«l.2Xt*C0EFRP:*.G9.‘»*2X* 

1*C0NSTA:*.G3.4.2X.*PRTCHG=*.39.4,2X. • GMASS = * . G9 .4 . 2X . • P RBDIA = » . G9. 

2*t ) 

WRITEt 6.7121 VELMAG.PBVELC.WEIGHT.VELVEC 

712 FORM AT (IX . *V EL KA G= * . G 9 .4 . 2X » *PBV ELC = * . G9 .4 1 2X • 

1*U EIGHT: *.G9. A. 2X * » V EL VE C= * . 3 ( G 9 . 4 . 2X )) 

CALL C0NTM2 5 REDO THE CC NT AM CALC. W/PRINTING. 

GO TO 600C 

36 CONTINUE 

00 37 1:1.3 

37 NEWVEL (I) : VELVEC(I) 

38 WNETt 4|:WNETt 4) ♦WEIGHT 

WNET (5>:WNET (5 I ♦SORT CN EW VE L 1 1 > ••2*NEWVELC2 ) * *2 +NE WV EL ( 3 ) * *2 ) 

1 *W EIGHT 

TPETA:-ATAN (NEWVEL (3J/SGRT (NEWVEL(1)**2'*NEWVEL(2 ) *«2*1 .E -9 1MPI/2 . 
TARG=THETA*17.39/PI*1. 
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WANGDS(IARG1=WANGDS*IARG!* WEIGHT 
GO TO 6000 

5999 CONTINUE 

6000 CONTINUE 8GRAIN VEL LOOP END 
5C0C CONTINUE 9GPAIN SIZE LCOP END 

IF( WN£T( 41.GT.0.1 WNET*5!=VNeTt5!/VNET*4! 

WRITEC6. 66661 
6666 FORMAT* IX* 15* /» IX 1 ! 

WRITE* 6*9100 IPOS HID 

9100 FORMAT* 1H0* *80X SCORE FOR RANDOM POSITION*. 2X .1P3E9.2 1 
WRITE * 6*9101 ! (WNET (I) *1:1*4) 

9101 FORMAT* IHO.’NO. SAFEHITS* *2X. G9.5*/*1X* * NO. RECONTAM HITS**2X*G9.5 
1»/»1X**N0. NEAR MISSES **2X»GS.5*/*1X**N0. ESCAPES** 2X»G9.5 ) 

IF*WNET* 9I.GT.0.I WPITE*6.9102» WNET (51 

9102 FORM AT * IX * • AVE • ESCAPE VEL. *»G9.5I 
00 92bO 1=1*9 

9200 WNET 2* I)=WNET2 *1! ♦WNET (I) 

WNET2* 5) =WNET2( 5) *WNET (51 *W NET 19) 

DO 9201 1=1*18 

9201 VANGOS* 11= VANGDSd! *WANGDS(I J/N2 
IFCSMTEST .LE.0.5) GO TO 9970 

40 CONTINUE 3 POSITION LOOP END 
GO TC 947S 

94’C WRITE* 6.94741 IMM . IN V « 3M TEST 

9474 FCRMAT<1HC»’~HE METECRCIC CLASS OF MASS C ROU P* » 16 . 2X * * AND VELOC. 
1CR0UP *. 16./ .IX. ’YIELDED C N LY* . G12 .4. 2X * » RELEAS ED GRAINS.* 

2./.lX.*GC CN TO THE NEXT MET ECROID CASE*! 

9475 IF* WNET2*4! .3T.0.1 VNET2 * 5) =W NEI2 «S> / WNET2 * 4 » 

WRITE* 6.930CI IHM.IMV.N2 

WR ITE t 6 » 6665 ! 

930r FORMAT * 1 H 1 . ’BOX SCORE FOR ONE METECROIC OF MASS GROU F * . 13 .2X » • AND 
1VEL0C. GRCUP*.I3»/»1X.* SUMMED OVER’. 16. 2X .• PCS IT IONS* I 
WRITE* 6.S1C1I *WNET2*I» ,I=X.4J 
IT* WNET2*4! .GT.0.1 WR ITE * o. 91021 WNET2*S) 

ACCUM1 = ACCUM1-*WNET2*1I*FLT0T 
ACCUM2=ACCUM2*WNET2*2)*FLTCT 
ACCUM3=ACCUMZ*WNET2<?!*FL’0T 
IVESC=IVESC*WNET2(41 *FLTOT 
SUMFL=SUMFL*FL T OT 

SMVESC=SMVESC*WNET2*5)*FLT0 T* WNET 2 *4! 

DC 94CC J=1 . 13 

9400 AN30IS* J )=ANGDIS*JJ *V ANGDS ( J ) * FLTOT 
WRITE I 6.9ECD I 
95CC FORMAT* 1 HI ) 

GO TO 5 

240 WR ITE* 6 » 241 ) 

241 FORMAT *1HC» ’FLUENCE TOO SMALL. SKIP THIS METEOROID*! 

GO TO 5 

5 F CONTINUE 

c • * • 

SAVESC=SAVESC/N2 
SAFEN0=ACCUM1/N2 
HITNC = ACCUM2/N2 
ANRMIS= A CCUM3/N 2 
IVESC= IVESC/N2 
SIGTOT=SORT*HITNO) 

WRITE*6.1CD! 

ICC FORMAT* 1H1.*B0X SCORE FOR OVERALL ENSEMBLE OF METEOROIO MASS AND V 
1EL0CITY GROUPS’! 

. WRITE*6. 94691 SUNFL 

9469 FORM AT t IX » *TCT AL NO. CF METEcROIOS IN ENSEMBLE DURING MISSION** 
1G10.41 

WRITE* 6*91C1 1 SAFENO.HITNO.AKRMIS.IVESC 
WRITE* 6*66651 
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WRITE<6*91C3> SIGTOT 

9103 FORMAT! 1H0. ’STANDARD DEVIATION FOR NO. OF REC OUT AMINAT I ON HITS IS*’ 
1.ZX.G9.5J 

RPR0B=1.-EXP!-HITN0) 

V AVESC=SMVESC/IVESC 
WRITEtS.760) VAVESC 

7GC FORM AT 1 1H0 * IX* ’THE AVERAGE ESCAPE VELOCITY IS • »2X . 1 PG9 .4 . 

12X * *M/SEC* ) 

DO 320 K=l»16 
320 ANGDEG! K>=10. 

WRITE (G»7G5) { ANGDEG! I I . ANGD IS (II >1=1*181 
765 FORMAT! 1H0* ’THE ANGULAR DISTRIBUTION FOR THE ESCAPED VELOCITY VECT • 
10R RELATIVE TO THE ♦ 7 D IRECT ION* .// * IX t »DEG . • *4X , »N0 . ESCAPES**//* 
?18t IX » CPF4 .0* 4X • 1PG10 • 3* / J J 

W D ITE I 6*335) RPROB 

335 Fo RMAT< 1H0 * 4 < IX */ ) * IX * *R ECO NT AKINAT ION PROBABILITY ISNGIZ.SJ 

c * * • 

STD? 

8Ct W^ITECEt eCl) 

901 FQRMATi IXt *THE MASS VALUE GMASS IS INCORRECT* I 
802 STOP 

81C UR ITE( 6 # 31 1 ) 

811 FORM AT C IX t *T HE VELOCITY VECT CR IS INCORRECT*) 

81? S T 0° 

END 
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c 


c 

c 


c 

c 

c 

c 

c 

c 

c 

c 

c 

c 


c 

c 


*•** ■« ADA JP- *•** 0. EDGARS (3 IONET ICS I « J. 3ARENG0LTZ (JPL) **•« 

REAL N'U. N £ X F 0 c * NEWVEL. IVFSC 

DIMENSION UNE7( 5) ,WN ET2( S> .VA n I (3E> , IS N(3C> *W A NODS (13 )» 

1 VANGDS ( 1 ” ) 

DIMENSION ANGDEG ( 1 9 ) 

DIMENSION PC~HIC ( 3 ) 

DIMENSION VA^dS). “ KACO C 3?) 

DIMENSION ANCOIF(ie). F LUX Ml C? » 1C 1 

DO M M CN 2LCKD1 CONNECTS CATAIN, MAIN AND RELEAS. 

CC“MON / ELC.X Cl / AKM.DPHC.K1.FR 2VEL (2E.1C) .CRflC )»SIG*DCIA (1C ) 

1. ANPART. ANN0RM.PR30IA (ID) 

COMMON 9LCKC2 CONNECTS CHKNIT. MAIN, RNDFCS AND SHADE. 

COMMON / 3L0XC 2/ AL(4) 

COMMON FLCKOI CONNECTS ESURF. MAIN AND RNDPCS. 

COMMON /3L0KC3/ ASP . ATP » ATS »= R I, AT OT . A3P , A3CS » All 

COMMON ELCKCA CONNECTS CCNTAH. MAIN, PARPOT, SPHIT AND TRAJEC. 

COMMON /3L0KC4/ COEFRP *CCNSTA , PRT CH3 , PHIO 
COMMON 3LCKCS CONNECTS DATAIN, MAIN, RELEAS AND YANC1 
COMMON /3LCKC5/ KC .®M ASStPR HO tR S ( 35) , FV AP (35 ) * FPKACC < 35 J . RHO 
COMMON BLCKCC CONNECTS DATAIN. MAIN AND YANG1. 

COMMON /3L0KC 5/ E.H.NU.PC.ISKIP 

common SLOKce connects ccntak. catain, main and parpct. 

COMMON /3L0KC8/ DIM »P ATHMX » SR » VELM A3 

COMMON 2LCKC? IN CHKHIT , CC NT AM , M AI N,S H A DE , S PH IT » T RA JEC AND TRANSL. 
COMMON /3LCKC9/ NEWP0S(3) 

COMMON DLCK1C CONNECTS CONTA H, MAIN. SFHIT .TRAJEC AND TRANSL. 

COMMON /3L0X10/ NEUVEl«3>. DELTAT 

COMMON SL0K11 IN CCN’’ AM . EFIE LD . MAI N. RNC PCS »S PHIT . TRAJEC .V NCRML. 

COMMON /3L0K11/ P0SVEC(3> 

COMMON PL0K12 CONNECTS CONTAK. MAIN AND TRAJEC. 

COMMON / 3L0K 1 2/ GMASS .VELVEC(3J »ILCST 

COMMON BL0K13 CONNECTS CHKHI 7 » CONTAK. MAIN AND SPHIT. 

COMMON /3L0K13/ IHITC2) 

COMMON BLCK1A CONNECTS CCNTAM. EFIELC, MAIN. SHADE AND TRAJEC. 

COMMON /3L0K14/ IP AN El . ISHADI 

COMMON BLCK1S I N CHKHIT .CO NT AM* EFIEtC. MAIN. RNDPCS. SHADE. VNORML. 

q.ommon /Blok is/ isectr 

COMMON SLCK16 CONNECTS MAIN. RNDPOS .SHADE AND VNORML. 

common /blokib/ npnl 

COMMON BL0K1 7 CONNECTS C A T A I N» MAIN AnD YAnGI. 

common /Blok it/ pvel.accmin 

COMMON BLOKie IN OAT AlN . E r IE LD *ESU RF . MAIN. RNDPOS .SC PLO T» SHADE . 

COMMON /SL0K1 3/ SP1 .SP 2 .SP 3 »SP 3 

COMMON BLOKlf CONNECTS DATAIN. EFIELC. ESURF. MAIN. RNDPOS AND SCPLCT. 
COMMON /3L0K13/ SP01.SPC2 .SP3 3 

COMMON BLCK2C IN CHKHIT. EFIELC. MAIN. RNDPOS. SHADE AND VNORML* 
COMMON /3L0K2C/ TAN T H1 »T ANTH2 »T ANTH3.T ANT H4 

COMMON EL0K21 IN CHKHIT. CO NT AM. EFIELC. MAIN. RNDPOS.SHAOE* SPHIT. SCPL 

COMMON / BLOK 21/ Zl. 22 . Z3. Z4 »Z 5 

COMMON BLCK22 CONNECTS MAIN AND VNORML. 

COMMON /BLOK 2 2/ VELNPM(3) 

COMMON BL0K23 CONNECTS EFIELC. MAIN AND SHADE. 

COMMON / BLOK 2 3/ XYUAKE.ZVI AKE » JPOT 

COMMON BL-CK2A CONNECTS EFIELC. ESURF AND MAIN. 

COMMON / 3L0K2 4/ AL AM3 12) » ANP HOT ( 2) . AKPHOT ( 2 ) . EFEL < 2 ). EFED ( 2 ) . OEB. 

1AM AT (9), EL (2 ) * PHIEL (2) .PHIED (2) . ALAMAV . PHIAVE 

COMMON '3L C< 2 2 CONNECTS CONTAM. EFIELO. MAIN AND TRAJEC. 

CCMMCN /DLC/.2S/ EVEC(I) 

COMMON CL OK 23 CONN EC T D DATAIN. MAIN AND FARPCT. 

COMMON /3LCK2C/ MAT * AV r 

COMMON 3LCX23 CONNECTS DAT A IN ♦ E c IELD , ESU R r « MAIN AND PARPOT. 

COMMON /SLCXD2/ AKT E . AKTP . ANE 

COMMON 3L0K3C CONNECTS CATAIN AND MAIN. 

CC M MON /PL0K3C/ACRH0 (?) . D1 . D 2 . C3 . Dt . Nl . N2 .S RMAX .S RMIN 
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1.EPSILN.ELEC7 » A" VEL 13) .TMISS. XPL OT 

: common sl'-kt: connects ccvam, efielc and main. 

COMMON /9LCX 31/ NRMIS5 

ECU IV AL CN CE (FVAP (1J .VAP (1) ) » ( F PKAC C (1) » PK ACC 1 1 ) ) 3 MAIN. RELE A s 

EQUIVALENCE (KlO.jSIZE) . C KcD. JVEL) 

: FORMAT I IHI •/// »1X » ’MASS GROUP ’ » 1 3 . 5 X . • V ELO C . GROU P • . 13 . 5 X • 
1’HISTORY NO.’ »I5 ) 

R EAC( S.16CC ) KML.KMU *MVL»MVU .IGML*IGMU .N’2. IHICUP 
WRITECo.IGCCIMML.MMU.MVL.MVU. 13 ML . 13 MU . N2 . IHI CU P 
1CCE FORMAT (8(13.3X1) 

00 355 IT-1 ♦ IHI CUP 

255 DUMMYrR ANCNC (1 ,.C.) S INCREMENT THE SYSTEM RANDOM NL'MSER GEN. 

Pl = 3. 14159 
DRHOiADPhO (MAT ) 

ALPHArATAN C2.*l.oE4*SQRTCAKT? )/ <AVP*S 3RT (Pin ) 3 WAKE HALF ANG. 

DRPIC6=CRHC*PI/6. 

TANTH1=D1/CS?4-SP03> 

TANTH2=-D2/ (5P4-SPD2) 

TANTH3=-C4/ CSP2-SP01 ) 

TANTH4=TAN (P 1/ 2 . - AL PH A ) 

700 F9RMATC IX. ’THE TANGENTS FOR CONE AN3LES ARE’ *4( G9.4.2X ) I 
WRIT EC 6.7CD) T ANT HI. T ANTH2.T ANTH3.T ANT H4 

0 

XYWAKE=CSP4-»SGRT CSP 4*SP 4*8 .* SP 1* CS P3-S P4 1 / PI 1 1/2 . 
ASP=16.*SP1*CSP3-SP2> 

ATF=PI*SPC1**2 

ATS=PI*C SPC1*SP2)*SQRT (04*04* (SP2-SP01) • CSP2-SP01H 
Al=ASP*ATS*ATP 

AI I=PI*( SP02*SP4)*S0RTCD2*D2* (SP4-SPG2) * (SP4-SP02I) 

ATOT=AI*AII 

PRI=AII/ATOT 

ABP=PI*SPC3**2 

ABCS=PI*t SP03*SP4)*S3RT(01*C1*(SP4-SP03) • CSP4-SP03I) 

WR iTe< 6.70 2) ASP . ATP . ATS . AI* A II* AT OT . PRlt ABP. ABCS 
702 FORMATC IX.* ASP = ’» G7.3.3X* »ATP =*»G7.3* 3X* • ATS = * * 07.3* 3X» « AI=«* 

1G7. 3.3X. *AII=*.G7.3.3X. *AT0T:*.G7 .3 . 3X » • PRI= • • G7 .3* 3X * » AB P= *. G7 .3 
2.3X.*ABCS=*.G7,3) 
sgtarg=epsiln* akte/ JANE*ELECT) 
sumfl - 0.0 

IVESC = O.C 
SMVESC=0. 

ACCUKlzC. 

ACCUM2Z0. 

ACCUM3zO. 3 FOR NRMISS. SEE EFIELO 

00 3 LM=1»1S 

3 angdisclm) =c . a initialize this array explicitly. 

10 SR=10.**<C.5*(ALCG10ISRMIN> ♦A L0G10 (S RMAXI ) ) 

ZlzSORTCSQT ARG)*SR a CEBYE LENGTH AT RAOIUS SR 

Z2=Z1*01 

Z3=Z2*D2 

Z4=Z3*D3 

Z5=Z4*D4 

D£8=Z1 

CALL ESURF 

WRITEC6.720)Z1.Z2.Z3.Z4.ZS 

72C Form at C IX. *Z1=*. G9.4.2X. *22 = S G9 .4. 2X. *Z3=*. G9 .4»2X» »Z4=*»G9.4» 
12X»«Z5=*.G9.4) 

ZWAKE=Z4-XYW AK E*T ANT H4 
WRITEC6.701) XYWAKE.ZWAKE 

701 FORMAT C IX . *X YW AK E= •* FS .4 . 4X. *2 W AKE= * . F9 .4 I 

C THE AL’S ARE THE Z-AXIS INTERCEPT FOR THE APPROPRIATE CONE. 

ALC l)zZ2-SP4*T ANTH1 
ALC 2)=Z2-SP4*TANTH2 
AL(3)=Z4-SP2*TANTH3 
AL C 4 )rZW AKE 
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WRITE* 6.721) AL 

721 FORMAT* IX* 'THE AL ARRAY CO NT A INS**4(G9.4.2X) I 

I FIKPLOT .EG .1) CALL SCF LOT * 3G OON) 8 SPACECRAFT PLOT. 

IF* IGOON.EQ.-l) STOP 

SM0LD1=1C.** (-18.142.9.33 - 1 .5* ALC61C *SR) ) 

DO 200 IMV=MVL»MVU 9 MVU.LE.3 

PVEL=APVEL *IMV J/SQRT (SR J 

DO 250 IMM=MML»MMU a MMU.LE.1C 

PMASSE =IMM-12 

PM0SSE=PMASSE-1. 

If'(IMM7e6.1> PMASSE=-10.52 
IF* IMM. EQ. 10) PMASSE=-2.52 

C CALCULATION CF COMETARY mETE^CIC FLUX PER NASA SP 8C83. 

IF(PHASSE.LE.-S.C) GO TO 11 
SLOGlO=-18.173-1.213*PMASSt-1.5* AL0G1 C (SrJ 
. GO TC 12 

11 SLC GlC =-18 .142-1 ,534*PM ASS E~ .0 63* FMASS E*PMASSE-1 .5* A LOGIC ( SR ) 

12 PMASS = 10.** (PMASSE- 3.) 3 IN KILOGRAMS. 

SMNEW1=1C.** (SL0G1C) 

IFCIMM.E3.1) GO TO 25 

IFCIMM.EQ.1CJ GO TO 25 

IFCIMM.EQ.2) PM0SSE--10.52 

IFCPM0SSE.LE.-G.) GO TO 13 

SL AG10:-13.17 3-1.213*PM0SSE-1 .5* A LOGIC (S R) 

5MCLDl = 10.** (SLAGIC) 

GO TO 25 

13 SLA GlC:- 18 .142-1 ,5 84*PMCSSE- ,C £3 • PM CSS E' PMCS S E-l .5 *A LOGIC (SR) 

SMOL 01=10.*** SL AGIO) 

25 SMlrSMNEWl 

I F( IMM .L E . 9 ) SM1 = SM0L01-SMNEW1 3 GROUP NO. DENSITY 

2E FLUXM1C IMV.IMM ) = .25* P V E L*S Ml /3 . 

FLT0T = FLUXM1( IM V. IMM) *AT0T*TMISS* 8.14 E4 
WRITE(6*7C5) IMM.IMV .PM ASS . FV EL* FLUX MI CIMV . IMM )» FLT CT 
70S FORMAT* IX. 'MASS GR0UP='.I3. 2X**VEL0C. GROU P=' . 13. /. IX. ' THE METEOR 
1CI0 MASS IS '.1PG1C.4.2X. 'THE METEOROID VELCC. IS'»1FG1C.4» 

2/. IX. 'WITH A GROUP FLUX 0 F ' , 1 ° GIG .4. 3X.'AND 

3.TFE TOTAL M ETECRC ID IMPACTS CF THIS TYPE EXPECTEO ARE’.G9.4) 

I^( FLTOT.LT. 0.01) GO TO 24C 

CALL YANC1 3 CESCPIEES SURFACE MV’MT 

IF ( ISKXP.E3.1 ) GO TO 75C 

CALL RELEAS 3 CET . F A RT . RE LE A SE PARAM. 

DO 80C1 1=1.18 
8CC1 V AN GDS ( I) = 0.0 
DO 8000 1=1.5 
8CCC WNET2C I)= 0.0 
C ALGORITHM TO PACK VA° ARRAY 

V AP 2 C 1 ) = ABSCVAPC1)) 

I J=2 

ISFUN ( 1 ) =1 
DO 8100 1=2. KO 

DC 82CC II=IJ.KC 
ISFUN* I)=II 

IFCABSCVAP *11) ).LE. 0 .5* ABS * V A P* IJ-1 III GO TO 83C0 
8200 CONTINUE 
83CC I FIX =ISFUN *1) 

VAP2* I)= 0.5* *ABS*VAP *IJ-1) ) .ABS*VAP* IFIX) ) ) 

IQU IT=I 

IF* ISFUN* I) .EQ. KO) GO To 8430 
IJ=lSFUN*l)*l 
810C CONTINUE 
840C CONTINUE 

IQUIT=IQUIT*1 
VAP2*IQUIT)=ABS*VAP*K0) ) 

WRITE* G. 3501) IQUIT 
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CO 35CO ITEMPrl.IQUIT 

WRITE* 6. 3502) IS FUN ( ITEMP ) » VAP2 *ITEMP> 

35CC CONTINUE 

3501 FORMAT* 1 HC . • I3U IT = » . 12 .// . 3X. * IS FU N* , 1CX . »V AP2 • . // ) 

35C2 F0RMATIEX.I2.1CX. 1PG9.4) 

C EXIT WITH GROUP POINTING ARRAY IN IS FUN* I) * I QUIT IS NO. OF GROUPS. 

C AND VAP2* I J IS GROUF VEL. ARRAY 

DO 40 J=1.N2 a INNER MONTE CARLO INDEX 

DC 4C0C 1=1.5 
400C WNET* I) =0.0 

DO 4CC1 1=1.18 
4001 WANGDS* I)=O.C 

WRITE*6»2> IMM.IKV.J 

CALL RNDPOS a DETERMINE IMPACT LOCATION 

DC 29 K = 1 » 3 

POSHID*K)=POSVEC*K) 

29 NEWPOS*K)=POSVECCK» 

CALL SHADE 
ISHHID=ISHACE 
IPHIO=IPANEL 
ISEHID=ISECTP 

WRITE* 6.709) POS V EC. IS HAD E. IPA NE L. NPNL. IS ECT R 

7C9 FORM AT* IX. »T HE RANDOM POSITION VECTOR IS*. 

1 3*G3.4.2X) ./.•ISHADE=».I2.2X. 

2 , IPANEL = , »I2.2X. •NPNL = »»I2.2X, •ISECTR=».I2 I 
CALL VNORML 3 NORMAL UNIT VECTOR 3P0SVEC 

WRITE* 6.710) VELNRH 

710 FORMAT* IX. »THE NORMAL VECTOR FOR THIS CASE IS •.3*69.4.2X1) 

I F* ISHADE .NE .0 ) GO TO 31 
CALL EFIELD 
P NIC =P HIEL CJPOT) 

GO TO 32 

31 PHICrPHIED *1 ) 

32 CONTINUE 

P HIDE =PHIC 

DO 5000 JS 12 E= IGML • IGMU 3 JSIZE.EQ.K1D GRAIN LOOP 
DIM=DDIACK1D)*1.CE-G 

CO EFRP=-7 .E-6*P I*. 25* * DIM/ SR) ** 2 3 RADIATION PRESSURE FORCE IN SUN 

C0NSTA=2 .*PI*EPSILN*DIM 3 CONVEX 'PHI TO -Q 

GMASS=DRPI0S**DIM**3) 

DO 6000 JV EL =1 . IQU IT 3JVEL EGUIV. KOC .GRAIN VEL LOOP 
C CALCULATE R EL AT IV E PROBABILITY OF JSIZE.JVEL PAIR 

IF* JVEL .EQ.l) GO TO 90C1 
IF* JVEL .EQ. IOUIT) GO TO 9303 
I FIX=ISFUN (JVEL) 

IFIX1=ISFUN*JVEL-1) 

PBVELC=PRBDIA*JSIZE)* (PRBVEHIFIX. JSIZE)-PRBVEL*IFIX1. JSIZE )) 

GO TO 9002 

SCO 1 P BV EL C=PRBDIA*JSIZE)*FRBVEL*1. JSIZE) 

GO TO 9002 

9003 FBVELCz PR BD IA IJS IZ E ) • *1 ,-FP BVEL * KC. JSIZE ) ) 

9002 CONTINUE 

WEIGHT rPBVELC* ANNORM STCTAL GRAINS THIS SIZE AND VEL REMOVED 8Y MT 
VELMAG= VAP2* JV EL ) 

DO 3E NN =1.3 S VELOC . VECTOR *INITIALLY) 

POSVECt NN) = ?OSHID*NN) 

35 V ELV EC* NN ) =V ELM AG* V ELNRM (N N) 3 EUILC INITIAL VELOCITY VICTOR. 

PHIO=PHIDE 
PRTCHG = CONST A*PHID 
IS HADE= ISHHID 
IFANEL=IPHID 
IS ECTR= ISEHID 

WRITE* 6.711) JSIZE. DIM. CO EFRP .CONST A. PRTCHG.GKASS 
l.PRBDIA* JSIZE) 
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711 FORM ATIlHC./// »1X. * jS IZ E = * » 1 2 . 2X . • C IK= * f G9 .4 . 2X . • C OEF RP= * . G9 .4 .2 X» 
l ’CONST A: *.0 3. 4. ?X » •<> RT CH3 - * . 3 3 . 4 , 2X . * 5 MAJS = * . G9 .4 . 2X . • PRBDIA = ’ « G9 . 

24) 

URITEI6.712) JVEl . V -LM A C.P C VC LC ,'wE I3HT . V ELv EC 

712 FCPM AT I IX, ’JVEL=*.I2.2X.’VELMAG=’» GS.4.2X. ’ PEV E LC = ’ . C9 .4 , 2X , 
1’U£IGHT=’.G9.4.?X.’VELVEC=’»3(o9.4,2X>> 

1FIWEIGHT .LE. .ECU GO TC SS9S 

ifi gm ass .le.o .> 20 to 3co a debugging 

V''AX=SCRTI-1C.*CC£FRP*PATHK!X/GMASS) 5 s times the minimum escape vel. 
IF| VMAX .UT.VELHAG.AND.NPNL.N: .-1) GOTO 36 
lF<V£LVEC< 3) .LE.r..ANC.NPN L .NE.-U CoToBIC S DEBUGGING 
CALL CONTAM 2 TEST CASE TRAJECTORY ANALYSIS 
IFI ILOST.EQ.l) gC t C 33 
I F| IHITI II .EO .11 UNCI tl)=WNE* (U ♦ WEIGHT 
I Ft IHIT f 2 ) .EC . 1 J WN ET I 2 J -W NE 1 1 21 *W EIGHT 
I PI NRM ISS.E3.1J WNETI3) = WNFT (3J*W£IGHT 
CC TO ercc 
3G CONTINUE 

WRHEI 6.736I VMAX 

7 35 FORMAT! 1H0. ’GRAIN V ElO CIT Y • .1 X. ’EXCEEDED ESCAPE 

1CRITER2A. VrL=**2X»lFG9.4) 

DO 37 1=1.3 

37 NEWVELIII: VELVECII) 

3? WJETI 4)=UNC T t 4) ♦WEIGHT 

WNET I E JrWNC’ (S I ♦SGR T It: EWVE L ( 1 1 •• 2*NtLV EL 1 2 I* *2 «NE U V EL I 3 > * »2 ) 

1 * W E : G H T 

’ FET A = - A’ AN If E WV EL I 3 ) /SORT ( LEW VEL 1 1 >«• 2 ♦NEW V E L 12 )*«2 ) ) +PI /2 . 
IAR6=THETA*17 .99/PIM. . 

W AN 60S I IAR G ) =W AN 60S 1 1 AR G > ♦ WEIGHT 
GO TO 6000 

5399 WR IT El E.5998 ) WEIGHT 

5998 FORMAT! IX. ’WEIGHT TOO SMALL!’ .39,5. •) .LOOP SKIPPED’ ) 

60CC CONTINUE 8 GRAIN VEL LOOP ENC 
5000 CONTINUE 3GRAIN SIZE LOOP END 

IFIWNETUI .GT.C.J WNETI5) =W MET 151/NNET I4> 

WRITE! 6. 91001 POSHID.IMN.IMV 

91CC FORMAT II HI. ’BOX SCORE FOR RANDOM POSIT ION’ .2X , 1P3E9 .2 . 2X . * AND METE 
lOROIO MASS GROUP’. 13. 2X . * METE OR OID VELOC. GROUP’. 13) 

WRITEI £.9101) IVNETII) .1=1.4 ) 

9101 FORMAT! 1H0. ’NO. S AFEHITS’ »2X» G9 .5 ./ ,1 X. » NO. RECONTAM HITS’, ?x» G9.5 
l./.lX.’NO, NEAR MISSES ’»2X.GS,5./.1 X.’no. ESCAPES ’. 2X.G9.E ) 

IF| WNETI 41 .GT.C. ) WR ITE I 6.91 02) WNETI5) 

9102 FORM AT I IX. ’AVE. ESCAPE VEL. ’.G3.5I 
00 9200 1=1.5 

92CC WNET2ID=WNET2 II) ♦W NET (I) 

DO 9201 1=1,19 

92-CI VANGDSIDr VANGDS II) .WANGDS 1 1 ) /n2 
4C CONTINUE 8 POSITION LOOP END 
WNET2I 5»=WNET2 I51/N2 
WRITEI 6,9300) IMM.IMV.N2 

930C FORM AT 1 1H1 » ’ BO X SCORE FOR METE CR OIC MASS GROUP*. 13. 2X» ’AND VELCC. 
1GR0UP’. 13. /.IX. ’SUMMED 0 VER » , 15 . 2X, • P OS IT IONS • ) 

WRITEI 6,9102) IWNET2II) ,1=1. 4) 

WRITEI 6. 9102) WNET2IE) 

ACCUM1 = ACCUM1+WNET2 ID.FlTCT 
ACCUM2=ACCUM2»WNET2I2) ’FLTOT 
ACCUM3=ACCUM3^WNET2I3)*FlT0T 
IVESC=IVESC*WNET2I4)*FLT0T 
SUMFL=SUMFL* FLTOT 
SMVESC=SMVESC^WNET2 15) .FLTOT 
DC 94CC J=1 , 18 

940C ANGDISl J ) = ANGDIS (J) +VANGDS I Jl.FLTOT 
GO TO 25C 
750 WRITEI6.751) 
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751 FCRMATUHC. 'ZERO IN AC Cl LEP All ON ARRAY. THEREFORE SKIP TO THE NEXT 
2 OUTER LOOP* > 

GO TO 2 5C 

40 URITE(6.241) 

41 FCRMATUHC. ’FLUENCE TOO SMALL* SKIP THIS LOOP’ J 
50 CONTINUE 3METE0R0I0 MASS LOO 3 E NO 
OC CONTINUE EMETECRCID VEL LCOF ENC 

50 CONTINUE 
C * * * 

SAFEN0 = ACCUM1/N3 
UTNC =ACCUM2/N2 
ANRMIS= ACCUM3/N 2 
I VES C - IVES C/N 2 
WRITE(S.ICO) 

ICC FORM A t C 1 Hi . * DC X SCOPE FCR OVERALL ELSFMPLE OF METECRCIC MASS AND V 
1EL0CITY GROUPS’) 

WRITE* 6.91C1 ) S AFENC .HITNC . A KP MIS. IVES C 
V A V ESC = SMVESC/ SU M FL 
WRITE* E . 7CC » V AVEE C 

76C FORMAT! 1 HC « IX .’THE AVERAGE ESCAPE V E L OCIT Y, 13 *. 2X . 1 PC 9 . 4 J 

00 320 K=l»18 
320 ANGDEGIK)=10.*K 

WRITE«6.765J CANGDEGCI) • ANGDIS 111.1=1*18) 

7S5 FORM ATC 1 HO. *T HE ANGULAR DISTRIBUTION FOR THE ESCAPED VELOCITY VECT 
10R RELATIVE TO THE *Z DIRECT I ON* .// .IX . *OEG . *.MX» »N0. ESCAPES*.//. 
21 8 ( 1X.F4»0.4X » IP 610. 3./) ) 

ANGTOT =C m 
DO 333 J = 1 * 18 

333 ANGTOT = ANGDIS (v'l ♦ANGTOT 

URITEI6. 335) ANGTOT. IVESC S DIAGNOSTIC PULL LATER 
335 FORM AT ( 1HC . ’ANGTOT = * * G1C .2 . 4 X » »IV ES C= » » G10 .2 I 

c * * • 

STOP 

800 HR ITE* 6 . 301 ) 

801 FORM ATC IX . ’THE MASS VALUE GMASS IS INCORRECT • I 

802 STOP 

810 WRITECG.811) 

811 FORMAT! IX. ‘THE VELOCITY VECTOR IS INCORRECT*) 

812 STOP 
END 
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SUBROUTINE PARPOT 

SUBROUTINE PARPOT 

CALCULATION CF the TIME DEPENDENT CHANGE IN POTENTIAL (PHICNG) 

AT A GIVEN POTENTIAL (PHIC) 

COMMON 3LCK04 CONNECTS CONT AM . MAIN. PARPOT. SPHIT ANO TRAJEC. 

COMMON /BL0K04/ CO EFRP .CONST A . PR TCH t . P HID 

COMMON BLOK 03 CONNECTS CONTAM. DAT AIN. MAIN ANO PARPOT. 

COMMON /BL0K08/ DIM .PATH MX »S R » VELMA C 

COMMON 8LCK14 CONNECTS CONTAM. EFIELD. MAIN. SHADE AND TRAJEC. 
COMMON /BL0K14/ IP ANEL . ISH AD E 

COMMON 3LCK26 CONNECTS DATAIN, MAIN AND PARPOT. 

COMMON /BLOK 26 / MAT.AVP 

COMMON BLOK 27 CONNECTS PARPOT AND TRAJEC 
COMMON /BLCK27/ PHICNG 

COMMON BLOK 28 CONNECTS DATAIN AND PARPOT. 

COMMON /BLCK28/ A ( 3 > .B ( 3 ) . M A TP OS ( 3 I .SOLSPC { 1 85 .2 I . Y IELD 1 15 .3 I »YIEL 
1DPC IS. 2) 

COMMON BL0K29 CONNECTS D AT AI N. EF IELC . ESURF .MAIN AND PARPOT. 

COMMON / 3L OK 2 9/ AK TE . AKTP . ANE 
EQUIVALENCE t DI AM .DIM ) . < RDIST »SR> 

4 FORMAT( 1HC. ’ERROR IN PARPOT. PHIO EXCEEDS RANGE AND E ®J ALS *. 3X» E10 
1.41 

5 FORMATI 1H0. ’ERROR IN PARPOT. YIELOP ARRAY. CHECK YIELD DATA’) 


The MODIFICATION NECESSARY Tc include the secondary electron CONTRIBUTION 
TO THE PARTICLE TRANSPORT ANALYSIS. REF. ’74 SEMI-ANNUAL REPORT ( JBB) 


PHIMIN=-3.*AKTE 

IF( ISHADE.GE.O.OR.PHIO.GT.PHIMIN) GO TO 15 
PHICNG=0. 

RETURN 

15 ONEALPrl.C 

IFtPHIO.GT.O.O) GO TO 20 
ALPHA=(ACMAT )-B(MAT)*PHIO)/( AKTE-PHIOI 
0NEALP=1. -ALPHA 
2C CONTINUE 

Cl=t ANE/4. ) »SQRT < AKTE) »G .7E5 

C2=( ANE/4.)*AVP 

C3=t AN£/1S.)*SQRT(AKTP)«1.£E4 

IFfISHADE) 15C0.10.110C 

CONTINUE 

IFRONT r C 

TOTYLO = 0.0 

M ATP = MATPOSIHAT) 

PHIUSE = 0.0 

IF CPHIC .GT. 0.0) PHIUSE = PHIC 

PHIBAS IS THE WORK FUNCTION 

PHIBAS = YIELDP(MATP.l) 

ENGBAS = PHIBAS ♦ PHIUSE 

PHIBAS IS THE PHOTON ENERGY AT THE BOTTOM OF THE FIRST NON-ZERO 

YIELD RANGE. I.E. WORK FUNCTION (APPROX.) 

ENG3AS IS THE MINIMUM PHOTON ENERGY THAT MAY CAUSE 
PHOTOEMISSION (APPROX.) 

IF (ENGBAS .LT. .001241) CC TO 30C 
IF (ENGBAS .LT. 23. G) GO TO 1'50 
GO TO 2CCC 
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c 

C THE ENERGY IN SOLSPC CLOSEST TO ENGBAS IS FCUND ANO THE RELATED 
C INTEGRAL OF IRRADIANCE STORED IN SPEC 1* IF A HIGHER RANGE OF 

C YIELDP WAS NEEDED v THE YIELD FROM THE PREVIOUS PARTIAL RANGE IS 

C STORED IN TOTYLD. 

C 

150 DO 200 J = l# 1 35# 1 

IF CS0LSPC(J#1) .LT. ENGBAS) GO TO 20C 
JS = J 

IF l(SOLSPC(J#l) - ENGBAS) .GT . ( ENGBAS - SOLSPC CJ-1#1 }} ) 

2 JS = J-l 

SPEC1 r SO LSPC t JS # 2 ) 

JT = JS ♦ 1 
MATP r HATP ♦ 1 

IF CIFRONT .EQ. 01 GO TO 500 

TOTYLD = TOTYLD ♦ f CSOLSFC 1 1 1 2 I - SPEC1 I • YIELD (MATP-2 # 

2 HATH / I YIELDP f MATP- 2* 2) • 1.6E-16I 

GO TO 500 

200 CONTINUE 

C NEW VALUES OF ENGBAS ARE CALCULATED BY STEPPING UP THE YIELDP RANGES 
C UNTIL A LARGE ENOUGH ENERGY IS FOUNC TC FIT IN THE SOLAR SPECTRUH 

C RANG E. 

C 

300 II = HATP ♦ 1 

DO 400 1= lit 15 • 1 

ENGBAS - YIELDP f I#1 ) ♦ PHIUSE 
IF (ENGBAS .LT. .0012411 GO TO 400 
HATP = I 
IFRONT = 1 
GO TO 150 

4C0 CONTINUE 

GO TO 2100 
C 

C THE YiaDS FOR EACH ENERGY RANGE OF YIELD ARE ADDED USING THE 
C IRRADIANCE VALUES OF THE SOLAR SPECTRUH UNTIL THOSE ENERGIES 
C EXCEEO THE SOLAR SPECTRUH RANGE. 

C 

500 DO 80C 1= N ATP t 16# 1 

ENGBAS = PHIUSE + YIELDP ( I #1 ) 

DO 600 J= JT# 185# 1 

IF C SOLSPC( J#1 1 .LT. ENGBAS) GO TO 600 
JS = J 

IF ( (SOLSPCI J#ll - ENGBAS) .GT. (ENGBAS - 

solspc< j-i#im js = j-i 

SPEC2 ~ SOLSPC (JS #2) 

JT = JS ♦ 1 
GO TO 700 
CONTINUE 
GO TC 850 

TOTYLD = TOTYLO ♦ (SPEC1 - SPEC2) • YIELD! 1-1# HAT) / 

( 1.6E-1B • YIELDP f 1-1 # 21) 

SPEC1 = SPEC2 
CONTINUE 
ENGTOP = 24.8 - PHIUSE 
C 

C HERE A FINAL YIELO IS ADDED FOR THE HIGHEST ENERGIES OF THE SOLAR 
C SPECTRUM. 

C 

TOTYLD = TOTYLO ♦ ( (SPEC1 - S OLS PC 1185 1 2 1 ) • 2.0 • YIELDI15# HAT)) 
2 / f t 27* 5 ♦ ENGTOP) • 1.6E-161 

GO TO 900 

850 TOTYLD = TOTYLO ♦ C(SP£C1 - S OLS PC Cl 85# 2) ) • Y IELDI 1-1# MAT I ) / 

2 IYIELOP (1-1.21 • 1.6E-16) 


600 

700 

2 

80C 
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90C TOTYLD = TOTYLD • 1CCCC.0 / 2.0 
C 

C THE VALUE OF PHICNG IS CACULATEC USING THE VALUE FOR THE TOTAL PHOTO- 
C ELECTRIC YIELD OBTAINED — TOTYLD. 

C 

IF (PHIO .GT. O.C> GO TO 925 
C PHIO.LE.O. 

PHICNG-2#88E-9*3.14159*DIAM*( - C1»EX P C PHIO/ AKT E » •0NEALP*C2 ♦ 
1T0TYLD>/»RDIST*RDIST) 

RETURN 

C PHIO.GT.D. 

925 PHICNG - 2.88E-9 • 3.14159 * DIAM « (-C1 * SGRTtl.O ♦ <2.0 • 

2 PHIO/ AKTE )»*0NEALP*C2 ♦TOTYLD)/ ( RDIST *ROIST ) 

RETURN 

C 

C PHICNG IS CALCULATED WITHOUT ANY PHOTOELECTRIC OR PROTON EFFECT 
C BECAUSE THE PARTICLE IS IN THE SHADE AND IN THE WAKE. 

C 

1500 IF C PHIO .GT. 0.0) GO TO 1600 

C PHIO.LE.O. 

P HICNG = 2. 8 8E-9 *3.14 159*01 AM* (-C1 )*EXP( PHIO/AKTE )*ONEALP/ 
2(R0IST*RDIST) 

RETURN 

C PHIO. GT .0. 

1600 PHICNG=2 ,8 8E-9*3.14159*DIAm* (-01 J*SQRT (1.4 (2 .*PHIC /AKTE I > 
2*0NEALP/(RDIST*RDIST) 
return 

C 

c phicng calculated for shade eut outsioe the wake. 

C REF. *74 SEMI-ANNUAL REPORT ( J38) 

C 

HOC IF ( PHIO . GT.C. ) GO TO 1200 

P HICNG=2.88E-9*3.141S9*DIAM* (-C1 *EX P< PHIO/AKTE )*0NEALF*C3*SQRTC 
ll.-t2.*PHIO/AKTP>))/tROIST*R3IST) 3 PHIO.LE.O. 

I F(P HICNG .LT .C.O.ANC.PHIO.LE .PHIMIN) PHICNG=0. 

RETURN 

120 0 P HICNG-2 . 88 E-9 *3. 14 159* DIAM* (-01 »SGRT (1.4(2 .*PHIC /AKTE 1 1 •ONEALP* 
1C3*EXP( -PHI0/AKTP))/(R0IST*R3ISTJ 3 PHIO.GT.O. 

RETURN 

200C WR ITEt 6.4) PHIO 
CO TO 9C0 
21CC URITE(S.S) 

STOP 

END 


ORIGINAL PAGE IS 
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• ••• NASA JPL *••* D. COGAPS (BIONETICS)* J. BARENGOLTZ (JPL) • *•• 
VERSION OF MAIN FOR PYRO CALCULATION IN RECONTAMINATION 
REAL NU* NEWPOS. NEWEL. IVESC 

DIMENSION WN ETC 5) .UNET2C5) *VAP2 (36) • ISFU N(36 t «W ANQOS (18 1* 

1 VANGDS CIS) 

DIMENSION AN8DEGC18 I * AN A ME (3) 

DIMENSION P0SHXD(3) 

DIMENSION' VAPI 35)* PKACCC35) 

DIMENSION ANGDIS (18 ) 

C COMMON 3L0K01 CONNECTS DATAIN* MAIN AND RELEAS. 

COMMON /BLCKC1/ AKM .OR HC * K1»FR BV EL (35. 1CI t CR(10 I.S IG*DOIA (1C I 
1 » ANPART * ANNORM* PRBDIA(IC) 

C COMMON BL0K02 CONNECTS CHKHIT* MAIN* RNDPOS AND SNAOE. 

* COMMON /BL0K02/ AL(R ) 

C COMMON BLOK03 CONNECTS ESURF* MAIN AND RNOPOS. 

COMMON /8L0KC3/ ASP » ATP .ATS .PRI.ATOT * ABP . ABCS * All 
C COMMON BLOKCR CONNECTS CONTAN. MAIN* PARPOT. SPHIT AND TRAJEC. 

COMMON /BLOKCR/ COEFRP .CONST A *PPT CM3 » PHIO 
C COMMON BL0KC5 CONNECTS CATAIN* MAIN* RELEAS AND YANG1 

COMMON /8LOKC5/ KO»PM ASS.PR HO »RS ( 35) * FV APC 35 ) *FPKACC ( 35 )* RHO 
C COMMON BL0K06 CONNECTS DATAIN* MAIN ANC YANG1. 

COMMON /BL0KO6/ E.H.NU.PO.ISKIP 
C COMMON BLOKOB CONNECTS CONTAfi* DATAIN* MAIN AND PARPOT. 

COMMON / BLOK 08/ DIM .PATHMX. SR .VELMAS 

C COMMON BLOKOB IN CHKHIT » CO NT AM*MAlN*SHAOE»SPNIT .TRAJEC AND TRANSL. 

COMMON /BLOKO 9/ NEUP0S(3) 

C COMMON BL0K1C CONNECTS C0NTAh*MA1N*SpHIT.TRAJEC AND TRANSL. 

COMMON /BL0K10/ NEUVEl(3). OELTAT 

C COMMON BL0K11 IN CONT AM.EFIElO.MAIN* RNCPOS*SPNIT»TRAJEC»VNORML. 

COMMON /BLOK 11/ POSVECOI 

C COMMON BLOK 12' CONNECTS CONT AN* MAIN AND TRAJEC. 

COMMON /BLOK 12/ GmASS *VELVEC(3) .ILOST 
C COMMON BL0K13 CONNECTS CHKHIT* CONTAM* MAIN AND SPHIT. 

.• COMMON /BL0K1 3/ IHIT(2) 

C COMMON BL0K1R CONNECTS CONTAr* EFIElD* MAIN* SHADE ANO TRAJEC. 

COMMON /BL0K1R/ IPAnEL »ISHA0E 

C COMMON BL0K15 IN CHKHIT* CO NT AM*EFIE l0» MAIN* RNOPOS .SHAOE. VNORML. 

common /blok is/ xsectr 

c COMMON BL0K16 CONNECTS MAIN* RNDPOS .SHADE ANO VNORML. 

COMMON /BLOK 16/ NPNL 

C COMMON BL0K17 CONNECTS OATAIK* MAIN AND YAN81. 

COMMON /BL0K17/ PVEL.ACCMIN 

C COMMON BL0K18 IN DAT AIN. EFIE LD.ESU RF.MAIN* RNOPOS .SCPLO 1* SHADE. 

COMMON /BLOK 1 8/ SP1*SP2*SP3*SPR 

C COMMON BL0K19 CONNECTS DATAIN* EFIELD* ESU RF.MAIN* RNDPOS AND SCPLCT. 
COMMON /BLOK 19/ SP01.SP02.SP0 3 

C COMMON BL0K2C IN CHKHIT* EFIELD* MAIN* RNOPOS* SHADE ANO VNORML. 

COMMON /BLOK 20/ TANTH1*TAMTH2*TANTH3*TAMTHR 
C COMMON BL0K21 IN CHKHIT* CO NT AH*£FXELC* MAIN* RNDPOS .SHADE* SPHIT* 3CPL 

COMMON /BLOK 21/ Z1*Z2*Z3*ZR«Z5 
C COMMON BL0K22 CONNECTS MAIN AND VNORML. 

COMMON /BL0K22/ VELNRM(3) 

C COMMON BL0K23 CONNECTS EFIEtC* MAIN AND SHADE. 

COMMON /BL0K23/ XYUAKE.ZWAKE • JPOT 
C COMMON BL0K2R CONNECTS EFIELD* ESURF ANC MAIN. 

COMMON /JL0K2R/ALAM8(2) * ANPH0T(2) .AKPHOT (2)»EFEL(2)*EFED(2>*DBB* 

1AMAT ( 9) »EL C2)*PHIEL (?) .PHIEC (2) » ALAMAV . PHI AVI 
C COMMON 3L0K25 CONNECTS CONTAM. EFIELD. MAIN AND TRAJEC. 

COMMON /BLOK 25 / EVEC(3) 

C COMMON 9L0K25 CONNECTS D A T A I N ♦ MAIN ANC PA«?0 T . 

COMMON /8L0K 26 / MAT.AVP 

C COMMON BL0K2S CONNECTS DA T AIN . E^IELO . ESU RF . MA I N AND PARPOT. 

COMMON /BLOK ?9 / AKT E . A KTP » ANE 
C COMMON 3LCK2E CONNECTS DATAIN AND MAIN. 
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CCMMCN /3LOK3C/AORHC (?) . D1 . C 2* 03 • D4.N1 . NZ *S RKAX.C RMIN 
l.EPSILN. ELECT ,A»VEL <31 .TMISS. XFLOT 
C COMMON BLCK31 CONNECTS CC NT AM * EFIELC AND MAIN. 

COMMON / BLOK Z 1/ SS1KS 

C COMMON SL0K32 CONNECTS FYPCM ANC RnCPCS. 

C0MM0N/3LCK32/ °RMIN t PRM A X 

COUIVALENCE CPVAP (1) tVAF (1) > . (FPKACC (1 ) . PKACC (I ) ) S MAIN* RELEAS 
E3UIVALENCE ( KID. JS 12 E) * < KOO* JVTL) 

FIr3.141S? 

DRHO=ADPHOtMAT) 

ALPHA=ATAN |2.*1.EE4*SC.RT<AkTF)/C AVP*SORTIPI >>) a WAKE HALF ANG. 
DRPIOo=ORHC* n I/S. 

WRIT E ( E * 1 ) 

1 FORMAT! 1 Hl» *T HE FOLLOWlNS ARE GENERAL SPACECRAFT MODEL DATA’) 
TANTHlzDl/ C2F4-SFC3) 

TANTH2=-D2/ CSP4-SP02) 

TANTH3=-D4/ CSP2-SPC1) 

TANTH4 = TAN'(PI/2 .-ALPHA) 

WRIT EC E.7CC) T ANT HI *T ANT H2 * T ANT HI »T ANT H4 

700 FORM ATC 1 HC » * T HE TANGENTS FOR CONE ANCLES ARE* . 4 t G9.4. 2X ) ) 

C 

XYWAKE=C SP4+SQRT (SP4*SP4*3.*SP1* CSP3-S P4) / PI ) )/2. 

ASP=16.*SP1* CSP3-SP2) 

ATP=PI*SPC1**2 

ATS=PI*CSP01*SP2)*SQRT C D4* D4 ♦( S P2-S PCI > * CSP2-SPC1 )) 

AI=ASP*ATS*ATP 

AII=PI*CSP02 *SP4)*S3RT (02*02 + 11 P4-SPC2)* CSP4-SP02 I) 

ATOT=AI*AII 

PRI=AII/ATOT 

ABP=PI*SPC3**2 

ABCS=P I*ISP03*SP4)*SGRT(D1*D1*15P4-SPC3) *CSP>-SP03> ) 

WRITEC B.7C2) ASP* ATP* ATS* Alt All. AT OT. P.RI* ABP. ABCS 
7C2 FORK ATC1H0* • ASP = • . G7 .3 . 3X • * ATP = ♦ *G7 .3* 3X» » ATS = * * G7 .3 * 3X» • AI = «» 
1G7. 3.3X* • AII= *» G7.3.3X* *ATOT=*. G7 .3. 3X . • PRI=* • G7 .3* 3X * • ABP=*. G7 .3 
2.2X* *ABCS=*.G7.3) 

SaTARG=EPSILN*AKTE/ CANE* ELECT) 

SUMFL : C.C 
IV ESC = 0.0 
SMVeSC=C. 

ACCUM1=0. 

ACCUM2 -C , 

ACCUM3=0. 3 FOR NRMlSS* SEE EFIELO 
DC 3 LMrl.ie 

3 ANGDISCLM)=0. a INITIALIZE this array explicitly. 

1C SR = 1C.**(0.5*(AL0G1C(SRMIN) + JL0G1C(SRMAX>) ) 

Z1=S3RTC S3TAR G) »SR 3 0E3YE LENGTH AT RADIUS SR 

22=21*01 
Z3=22*D2 
24=Z3*D3 
Z5=Z4*D4 
DEB =Z1 

WRITE(G*720)Z1*Z2*Z3»Z4*Z5 

72C FORM ATI1H0. *Z1 = ». G9 .4 * 2X . »Z2 :* *G9 .4 . 2X * *Z 3 = • » GS .4 *2X . »Z4 = »*G9.4» 
12X ♦ *Z5= * * G9.4 ) 

ZWAKE = Z4-XYW AK E*T ANTH4 
WRITEC G* 7C1 ) XY WAKE* ZWAKE 

701 FORMAT C 1HC » ’X YU AK E = *» FS .4 *4 X. »ZW AKE=* » FS .4 I 

C THE AL * S ARE THE Z-AXIS INTERCEPT FOR THE APPROPRIATE CONE. 

AL ( 1 ) =Z2-S R4*T ANTH1 
ALt 2)=Z2-SP4*TANTH2 
ALC 3)=Z4-SP2*T ANTH3 
AL ( 4 ) =ZW AKE 
WR ITE C 6 *721 ) AL 

721 FORMAT! 1 HO* *T HE AL ARRAY CO NT A INS • , 4 C G9 . 4 * 2X ) I 
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CALL ESURF 

4 READ! 5* 1111 ) N2.KPL0T 
1111 FORMAT (2110) 

5 REA0I5»16C0» ERR =50 1 EN 0=50 ) HICUP.KO. TGML.I8MU. IPYRO. NPYRO. 

1PRMIN ,PRM AX • AN AME 

SMTEST=C. 

WRITEI6»1ECS) IPYRO 

1609 FORMATl 1H1. ’THE INPUT PYRO CASE DATA CARO FOR PYRO NUMBER*, 141 
WRIT El 6,1610 UHICUP, KOtIGML* I6MU,IPY RO, NPY ROtPRMIN, PRMAX , AN AME 
1600 FORMAT! 61 13, 2X) ,2E12.5,3A6) 

161C FORM AT 1 1HC. 6 113,2X1 ,2E12.5,3A6» 

IF<KPL0T.EQ.1 ) CALL SCPLOTIIGOON) 3 SPACECRAFT PLOT. 

lFlIC0CN.E3.-l) stop 

FLT&T=NP"Yr6 

WRITE(6,7C5) IP YRO » AN AM E»NPYfiC 

705 FORMAT! 1 HI, *P YRO NO . * »1X , 13 ,2 X, ’ LOCAT EO IN FIRST X~Y HjADRANT AT*, 

1//»2 7X,3AG»//»EX, ’TOTAL NO. AT ESUIVALENT POSITIONS ON S/C TS*»I4I 
WRITE! 6, 1593) 

5QC DC 51C I=1»KC 

READ! 5,16C1 ) IDUMYitRS (I) ,P KACC ( I) .VAPII) 

WRiTeI6,1601)IDUMY1,RS (I),FkACC 1I) tVAFII) 

510 CONTINUE 

1599 fORMaTI1hC.////»2X» *1 • , 11X , *R ! I) «, 9X , • PKACC ! I )*» 9X , •VARII )• ) 

16C1 FORMAT! I3.3X, 1P3E15.E) 

27 CALL RELEAS S DET .PART .RELEASE PARAM. 

DO 3001 1=1,13 
8001 V AN GDS 1 1 ) r 0 .0 
00 80CC 1=1,5 
8C0C WNET2I I)= O.C 
C ALGORITHM TO PACK VA’ ARRAY 

VAP2I1)= A3S !VA P (1J ) 

IJ = 2 

ISFUNI 1) =1 
00 8100 :=2,KC 

DO 8200 IlrIJ.KC 
ISFUNI I ) =11 

IF! ABS IV AP III) ).LE. C .5* ABS I V A FI I J- 1 ) ) ) GO TO 8300 
82CC CONTINUE 
83C0 IFIX=ISFUN II) 

V-AP2I I)= C.5* I A 3 S I V AP 1 1 J- 1 ) I+A3SIVAPI IF IX) ) ) 

I CU I T= I 

TFI ISFUNI I) .£2. KO) CO TC 843 C 

I J=ISFUN 1 1 ) •* 1 
8100 CONTINUE 
8400 CONTINUE 

IQUIT=IQUIT*-1 
V AP2I IaUlT ) = ABS IVAP I KO ) ) 

EXIT WITH GROUP POINTING ARRAY IN ISFUNI I) ,I3UIT IS NO. OF GROUPS* 

AND VAP2II) IS GROUP VEL. ARRAY 
DO 955 IT-1* IHICUP 

955 DUMMY=RANDNO (l.,0.) 9 INCREMENT THE SYSTEM RANDOM NUMBER GEN. 

00 40 J=1»N2 9 INNER MONTE CARLO INDEX 

WRITEI6,2) IPYRO, J 

2 FORMAT! 1H1,///,1X,*PYR0 NO . ** 13 , 5X» • P OS IT ION NO.*, 15) 

DO 4000 1=1,5 

4000 WNET 1 1) =0.0 

DO 4CC1 1=1,18 

4001 WANGDS! I)=0.0 

CALL RNDPCS S DETERMINE IMPACT LOCATION 

DO 29 K=l, 3 
PCSHIDlK)=POSVEC |K ) 

29 NEWPOSlK)=POSVEClK) 

CALL SHADE 
ISHHID=ISHADE 
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IPHID=IPANEL 
IS EHID= ISECTP 

CALL VNORKL a NORMAL UNIT VECTOR 3P0SVEC 

WRITEt 3.7C3) PO S V EC. ISHAD E. IPA NS L» NPNL. IS ECT R 

7C9 FORMAT 1 1HC • ’THE RANDOM POSITION VECTOR IS** 

131 G9. 4.2X) .// *1X* • ISHADE=».I2 *2X. 

2*IPANEL:**I2.2X. *NF NL = * * 22 * 2 X » *IS ECT Rr • * 12 ) 

WRITEt E * 71C ) VELNRM 

71C FORMAT tlHC . *THE NORMAL VECTOR FOR THIS CASE IS * . 3 (C1C .4 *2X )l 
IF< ISHADE.NE.O) GC TO 31 
CALL EFIELD 
PHIC = PHIEL< JPOT ) 

GO TO 32 

31 PHIO=PHIEDtl> 

32 CONTINUE 
PHIDE=PHIC 

DO 5CCC JSIZE=IGML *IGMU 3 JSIZE.EO.K1D GRAIN LOOP 
DIM=DDIA(K1D) *1 .OE-G 

C0EFRF=-7.E-6*PI*.2S* <DIM/SR)*»2 3 RADIATION PRESSURE FORCE IN SUN 

CONST A=2. »P I*EPSILN*D:M 3 CONVERT -PHI TO -3 

GMASSrDRPIC5»(CIM««3) 

DO SOCO JVEL = 1* I3UI"’ 3 JVEL F3UIV. KCD .GRAIN VEL LOOP 
C CALCULATE RELATIVE PRC3A3ILITY OF JSIZE*JVEL PAIR 

TFfJVEL .EO.i) GC TO 9001 
iFtJVEL .EQ. IOUIT) GO TO 9CC3 
IFIX = ISFUNM JVEL ) 
iFlXlriCFUN (JVEL-1) 

P3VEL CrPRBDIA tj SIZE) * IPR3VELI IF IX * JO I7E) - PR3V EL l IF IX 1 1 JSIZ E ») 

GO TO SCC2 

3Cri P3VELC=PRDD TA<J SIZE) *PR3 VEL (1 .JSI2E) 

cc to tree 

9C0 3 <>3vElC= P:-3DIA(JSIZE) * ( 1 .-? R 3 VE L( KO . Js IZE ) J 
9002 CONTINUE 

WEIGHT:? 2 VEL C*ANNCRM 3TCTA L 0 RAINS ■ T HIS SIZE AND VEL REMOVED 3Y MT 

vn m ac : vaf: <jvel> 

DO 35 NN=1.3 3 VELOC. VECTOR ( INITIALLY} 

FOSVEClNN)=PCSHICtNN> 

35 VELVECI NN ) iVELM AG* VELNRM (NN) 3 3UILD INITIAL VELOCITY VECTOR. 

P HD =P HIDE 

PRTCHG=CONSTA*PHIC 

ISHADE :ISHHID 

IPANEL=IPHID 

ISECTR =ISEHID 

AT EST: WEIGHT * FL TOT 

SMTEST=SMTEST-*ATEST 

I Ft WEIGHT .LE. .001. OR .ATEST. LE.0.05) GO TO 5999 
IF( 3MASS.LE.C.) GC TO 8CC 3 DEBUGGING 

VMAX=S0RT(-10 .*COEFRP*P ATHMX/GMASS) 3 5 TIMES THE MINIMUM ESCAPE VEL. 
IF(VMAX.LT.VELMAG.ANC.NPNL.NE.-l) GC TO 36 
IF( VELVECt 3) . LE.O.. AND. NPNL. NE.-1J GC TO 310 3 DEBUGGING 

CALL CCNTAM 8 TEST CASE TRAJECTORY ANALYSIS 
IFt NRMISS.E3.1) WNETt 3) = WN£ T 13) * WEIGHT 
IF( ILOST.Efi.il GO TO 38 
IF( IHITt D.EG.l) WN ET 1 1 ) :UNET tl ) + WEIGHT 
IFt IHITt 2) .EO.I) WNET t 2 ) :W NE T t 2 ) ♦WEIGHT 
IFt IHITt 2). EG. 11 GO TO 3CCC 
GO TO ECOC 
3000 CONTINUE 

c this prints the trajectory of the frcueeding recontamination event. 

DO 3500 NN=1»3 
POSVECtNN)=PCSHIDtNN) 

3500 VELVECt NN)-VELMAG*VELNRMtNN> 3 3UIL0 INIT IAL V ELOCITy VECTOR. 

phio=fhide 
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PRTCHG=CONSTA*PHIC 

ishade=ishhic 

ipanel=iphid 

tsectr=isehid 

WRITE* 6.2CCO) J. JSIZE. JVEL 

2 CCC FCRMAT(lHl»’ c OSlTlON*. 16. 4X. ’GRAIN SIZE GR0UP’»I6.4X. 

1 • 3 RAIN VELOC. GROUP’. 161 
WRITE* 6.711) DIP. CCEFRP .C ONST A . pRTCHG. GMASS 

1.PR30IA* JSIZE) 

711 FORMAT* 1 hC»/// .IX. *DlM=*.G9.4.2X»*C0ErRp=*.G9.4»2x» 

1’ CONST A=’. G3.4. 2X . »PRT CHG= ’ » 3 9. 4. 2X. » GMASS =* . 89 .4. 2X * * PRBDIA:*, G9. 
24) 

WRITE* 6.712) VELMAG.PBVELC.WEI3HT.VELVEC 

712 fO^MaT*1X* ’VElMAG=’.G5.4.2X»*pBVElC=*»G9.4»2X. 

l’WEIGHT=*»G9.4» 2X . • VEL V£C= ’ . 3 *G9.4.2X)) 

CALL CCNTM2 S REDO THE COnTAM CALC. W/pRINTlNG. 

GO TO SOOC 

36 CONTINUE 

DO 37 1=1.3 

37 NEWVEL * II = VELVECtI) 

38 WNET*4)=WNET*4) ♦WEIGHT 

WNET ( 5 l=WN ET *5 ) ♦SORT <N EW VEL * 1) ** 2 + NEWV EL* 2 )«*2«NEWVEL*3 )**2 I 
1»W EIGHT 

THETA=- AT AN IKE WV EL * 3 ) /SGRT * NEWVEL * 1 ) ** 2-»NEW V EL *2 >**2-»l .E-9 ) )-»PI/2. 
IARG=THETA*17 •99/PI+l* 

W ANGDS ( IARG ) =W ANGDS JIARO+ WEIGHT 
GO TO 6000 

5999 CONTINUE 

6000 CONTINUE SGRAIN VEL LOOP END 

5 CCC CONTINUE 6GRAIN SIZE LOOP C KC 

IF* WNET* 4I.3T.0.) WN ET * 5 ) =WNFT < 5 ) /W NE T * 4) 

WRITE* 6.6666) 

6666 FORMAT* IX. 15*/. IX) ) 

WRITE* E.910OPCSHID 

3100 FORMAT* 1HC. ’BOX SCORE FOR RANDOM PCS IT I ON’ . 2X . 1 P3E9.2 ) 

WRITE* 6.9101) (WNET *1) .1=1.4 ) 

9101 FORMAT* 1H0. ’NO. S AFEHITS ’ . 2X . 33 . 5 ./ . 1 X . * NO. REC ONT AM HITS’ . ?X . G9 .5 
1./.1X.’NC. NEAR MISSES ’.2X.CS. 5. /.IX. ’NO. ES C A PES • » ZX . G9 . £ I 

IF* WNET* 4I.ST.0.) WR ITE * 6. 91 02 ) WNETIS) 

9102 FORM AT * IX . ’ AVE • ESCAFE VEL. ’.G9.5) 

DO 9200 1=1.5 

92CC WNET2«II=WNET2 (I)*WNET(I) 

DO 9201 1=1.13 

9201 V AN CDS * I> - V AN COS C I> +W ANGDS * I ) / N2 
IF*SMTEST.UE.0.5) GO TO 9470 
4C CONTINUE 6 POSITION LOOP EnC 
GO TO 9475 

947C WRlTEt6»9474> IPYRO.SMTEST 

9474 FORMAT* 1H0. ’THE PYRO N DM 3ER ’ . 16 . 2 X ./ * IX . • Y IELDED 0NLY’.G12.4. 

izx. ’R eleased grains. ’./.ix.’cc on to the next pyro case*) 

9475 WNET2* 5) =WN ET 2* 5)/N2 

WR JTE * 6 ♦ S3CC ) IPYR0.N2 
WR ITE* 6* S66S) 

93CC FCRMATdHl, ’BOX SCORE FCR ONE FYRO OF TYPE N0.’»I3.2X. 

1/ . IX . ' SUMMED OV ER’» 16. 2X. ’POSITIONS’) 

WRITE* 6 . 9IC1 ) (WNET2II) .1 = 1.4) 

WRITE*6» 9102) WNET2 * 5 ) 

ACCUM1 = ACCUM1*WNET2 *1 )*FLTCT 
AC CUM2 = A CCUM2 ♦WNET2*2)*FLTCT 
ACCUM3 = ACCUM3+WNET2 * 3)*FlTCT 
IVESC=IVESC+WN£T2C4)*FLT0T 
SUMFL=SUMFL+FLTOT 
SMVESC=SMVESC+WNET2(5)*Fi_T0T 
DO 34CC J=1 . 12 
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34CC ANGDIS) J >=ANGDIS)J) *VAN GDS ( J > • FLT OT 
WRITE ( E.95CC ) 

9500 FORMAT ( 1 HI ) 

GO TO E 
5C CONTINUE 
C * * * 

SAFEN0=ACCUM1/N2 
HITNO = ACCUM2/N2 
ANRMIS= ACCUM3/N 2 
IVESC= IVESC/N2 
SIGT0T=S3RTIHITN0> 

WRiTECGfICOJ 

IOC FORMATi 1H1» 'BOX SCORE FOR OVERALL ENSEMBLE OF PYROS*) 

WRITE ) E » 946 S ) SUMFL 

9<*G9 FORMAT) IX. 'TOTAL NO. OF PYROS IN ENSEMBLE OURING MIS SION** 610.9 ) 
WRITECE»91CU SAFENO .HITNO » A NR MIS* IV ESC 
WRITE( E* 6665) 

WRITE) 6.9103) SIGTOT 

3103 FORM ATI 1 HO. 'STANDARD DEVIATION FOR NO. OF RECONTAMINATION HITS IS* 
1.2X.G9.E) 

RPR03=i. -EX?) -HITNO) 

VAVEEC r SMVESC/ SUMFL 
WR ITE) 6 » 7SC I VAVESC 

7EC FCRM AT { 1 HC * 1 X * 'THE AVERAGE El C A PE VELOCITY IS ' . 2X » 1 P09 .9 . 
12X.'M/SEC*> 

t DC 320 K=1»1E 

320 ANGDE G) K ):10. *< 

WRITE«6*76E) ( AN GDE S( I) » AN GDIS (11*1:1*18) 

765 FORMAT) 1HC* 'THE ANGULAR DISTRIBUTION FOR THE ESCAPED VELOCITY VECT 
10 R RELATIVE TO THE *7 D IRE CT I ON' ./ / , IX . 'DE G . ' . 4X * ♦ NC . ESCAPES*. //, 
2131 IX .CPF4.0. 4X . 1PG1D .3./J ) 

WRITE ( E.32E1 9 FRO B 

335 FORMAT) 1 HC . 4 l IX . / ) . IX. 'R E CO NT A "INAT I ON PROBABILITY IS'.G12.5) 

c * * t 

STOP 

8CC WRlTE(6»eri) 

801 FORM AT< 1 X • # THE MASS VALUE GMASS IS INCORRECT* ) 

802 STOP 

810 WR ITE( £ ♦ 811 ) 

811 FORM AT ( IX t f T HO VELOCITY VECT CR IS INCORRECT*) 

812 STOP 
END 
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function randnocbtmlim»uprlik> 

• *** NASA JPl 9/12/74 *••* 0. EDGARS* BIONETICS ,*••• 

**** 

RETURNS A RANDOM NUM3ER UNIFORMLY SELECTED FROM THE INTERVAL 
(STMLIM.UPRLIM). DUMMY IE A SUPERFLUOUS VARIABLE 

**•* 

R=UNIFRMCDUMMY ) 

RANDNO=BTMLIM*U.-R)*UPRLIM*R 

return 

END 
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SUBROUTINE RELEAS 


SUBROUTINE RELEAS 

• *** NASA UP L 9/12/7«* ***• C. EDGARS. EIONETICS • ••• 

THIS SUBROUTINE IS THE AO APT I ON OF THE PARTICLE RELEASE 3TU0Y SOFTWARE 
DEVELOPED AT JPL BY C. BAUERLE AND J. B. BARENG0LT2 11973). 

THE OUTPUT IS IN THE FORM OF TWO THRFF-OINENS I0NAL ARRAYS OF DATA 
corresponding to the hctc NM tv input parameters. 

* * * • 

INPUT PAE AMEIERSJ _CFCR SURFACE BEHAVIOR AFTER METERCIO IMPACT) 

1 ) RS - AN ARRAY OF RADIAL DISTANCES FOR THE VELOC. AND ACCEL. ARRAYS 

2 ) PK ACC - THE PEAK ACCEL. S A DISTANCE RS FROM IMPACT 

3) VAP - THE SURFACE VELOC. 9 ACCEL. PEAK. 

A ) KO - THE RS* PKACC* AND VAP ARRAY DIMENSION. 

OUTPUT parameters: , 

1 ) DDI A - THE ARRAY CF EJECTA DIAMETERS ISC-1CC MICRONS). 

2) K1 - THE ARRAY DIMENSION OF OOIA. 

3) CPIK1) ARE THE CLEARING RADII. 

4 ) PR8VEL • THE ARRAY OF EJECTION PRC3A3ILIT IES .VS . 1 VAP.ODIA ) VALUES. 


DIMENSION VAP I 35 ) 

DIMENSION R(35 )*A(35)*RF(35*10)*AK (1C. 35) .PVEL 135*10 
COMMON BL0X01 CONNECTS DATAIN* MAIN AND RELEAS. 

CCMNON /BLOKOl/ AKH *DRHC *K 1* FR 8VEL 1 35. 10) . CR (1C ) .SI6.DCIA <1C ) 
1* ANPART • ANNORM.RRBOIAllC ) 

COMMON BL0K05 CONNECTS DATAIN* MAIN* RELEAS AND YANG1 
COMMON /3L0K05/ KO.PM ASS.PRHO .RS135J .FVAP135) .FPKACC1 35 ) . RHO 
EQUIVALENCE (PRBVELl 1*1 ) »PVEL Cl *1 1 ) * (RS 1 1) *R 1 1 ) ) 

EQUIVALENCE 1 FVAP (1) * VAP 1 1) ) * (FPKACC (1) . All ) ) * 1 INAXI. KO ) 

NUMBER OF CASES (M) MUST BE .LE. NH 

DEFINE DIS(D)=1./(D*0) 

ANNORM: 0.0 
IMAXII=IMAXI-1 
PI=3. 14159 

Cl=(PI*0RH0*9800»*01.E-12)/6. 

S3RT2=S0RT(2.) 

CALCULATION OF PERCENT REMOVED 

00 1 11=1*10 
DDIAl II)=II*10 
DD=DDIA(II)*DDIA(II)*C1 
00 2 IJ=1* INAXI 
Y=DD*A (IJ) 

AK1 II* IJ )= AL00101Y) 

X=ABS( AKM-AK (II»IJ))/tSQRT2*SI6) 

IFIAK(II.IJ) .LT. AKN) 00 TO 3 
RF1IJ.II) =0.5*11. ♦ERF(X))*1CC. 

00 TO 2 

3 RFf IJ »II) = 0.5*ERFC(X)*100. 

2 CONTINUE 
1 CONTINUE 

CALCULATION OF CLEARING RADXUS.CRtH.) AND RELEASE vriociTY 
PROBABILITY* PIV .OT. VO) 

DO 5 11=1.10 
CR1 II) =Pt 1 )**2 • 

PVELl l.II)=CR(II) 

DO t I J=1 * IM AX II 

CR1 II) =CR( II) *RF 1 I J* II) *R (IJ) *(R1 IJ*l)-Rl IJ) )/50. 

G PVEL (IJ*1*II)=CR (II) 

CR(II)=SQRTICR(II) ) 

DO 8 K=1»IMAXI 

8 PV ELI K • II) =PV EL IK * II) / (CR (II) *CR(II) ) 
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E CONTINUE 

00 1CCC 1=1*10 

rOBOlA(I)=Pl*CRlU»CO(I)*CCIS(DDIACIl-E.J-DIStDOIA(I)'*5»>) 

i*ANPART/0IS<5.» 

Ar;NCRM = ANf.'ORH + PRBCIA(I) 
lCCr CONTINUE 

DO 2CCC I=1»1C 

2000 PR3DIA! I ) =PRBOI A ( I ) / ANN9RM 
RETURN 
END 
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SUBROUTINE RNDPOS 


SUBROUTINE RNDPOS 

NASA JPL 9/12/74 ••** D. EDGARS* BIONETICS 

THIS ROUTINE PRODUCES A RANDOM POSITION ON THE SPACECRAFT AS DEFINED. 

X * Y * Z ARE THE POSITICN COORDINATES. 

IS IS INDEX IDENTIFINO SEC T SR OF IMPACT (1 OR 21 

NP IS 1 FOR TOP PANELS AND -I FOR THE UNDERSIDE OF THE SOLAR PANEL 
• • • • * • 

COMMON BL0KC2 CONNECTS CHKHIT * MAIN. RNCPOS AND SHAOE. 

COMMON /BLOKC2/ AL C 4 1 
C COMMON BL0K03 CONNECTS ESURF . MAIN AND RNDPOS. 

COMMON / BLOK 0 3/ ASP » ATP . A TS R I . AT OT . ABP . ABCS • All 
C COMMON BLCK11 IN CONT AM . EFIE LD . M AI N . RND POS . S PH IT *T R A JE C * V NO RML . 

COMMON / BL OK 1 1/ P0SVECC3) 

C COMMON BL0K16 CONNECTS CO NT AM.MAIN.RNCPOS .SHADE AND VNORML. 

COMMON /3L0K1 6/ NPNL 

C COMMON BLCK18 IN D AT A I N . EF IE LC *E SU RF .M A IN . RNDPOS . SCPLO 1. SH AOE . 

COMMON / 9LOK1 9/ SP1.SP2.SP3.SP4 

C COMMON BLOK19 CONNEC T S D AT AI N. EFIELC. ESURF. MAIN. RNDPOS AND SCPLOT. 

COMMON / BLOK 1 9/ SPOl » SPC2 . SPD 3 

C COMMON BLOK2C IN CHKHIT. EFIELD. MAIN. RNDPOS. SHADE AND VNORML. 

COMMON /3L0K2C/ T AN THJ , T ANTH2 . T ANTH3 , T A NT H4 
C COMMON EL0K21 IN CHKHIT .CO NT AM. EFIELC.MAIN.RNC POS. SHADE. SPHIT.SCPL 

COMMON / BLOK 2 1/ Z 1 . 72 . Z3 . Z4 . Z 5 
EQUIVALENCE (NP NL »NP 1 . ( IS ECTR .1 S 1 

EOUIVALENCE(POSVECd) . X ) . (POS VE C ( 2 ) . Y 1 . ( POS VEC ( 3 ) . Z 1 
PRrRANDNO(D. 0.1,1 
WRITE! 6. ini )PR 
1C1 FORMAT (1HC. , PR=*.G9.41 
NP=0 

IF(PR.GE.PRI1 GO TO 1C ii ( DETERMINE SECTOR FOR EVENT 1 

GO TO GO 

ID CONTINUE 5 IN SUN CR UNDER PANEL (WEIGHTED BY THE EXPOSED SURFACE! 

I F ( PR. L E . ( PRI *ATP / ATO T ) ) SO TO 20 

IF(PR.LE.(PRI.(ATP*ATS)/ATCT 11 GC TO 3C 

IF( PR.LE.(PRIM ATP»ATS*ASP/4. 1/ATOT) ) GO TO 40 
IF(PR.LE.(PR:-*(ATP*ATS-.ASF/2 .l/ATCTl! GO TO 45 
IF(PR.LE.(PRI*( ATP+ATS.3.*AS?/4.1/AT0T1 1 SO TO 50 
IF(FR.LE,(FRI.(ATP4ATS.ASP1/AT0T 11 GO TO 55 

GO TO 70 

2D X;RANDNC(r..SPCl! S TCP FLAT SECTION 

AY=SQRT( SP01*SP01-X*X ) 

Y=RANDN0(0.. AY1 
Z = Z5 
NF = 1 

GO TO 70 

3C ETA=RANDNO (C ,. 1.571 S SECTOR 1 CONE AREA 

RAD=RANDN0(S?01.SP2> 

X=RAD*C0S(ETA1 
Y=RAD*SIN( ETA ) 

Z = SSRT ( X * X + Y* Y1 »T ANTH3 ♦ AL < 2 1 

GO TO 70 

4 C X=R ANDNC(SP2.SP31 5 AEOVE SOLAR PANEL (11 

Y=RANDNO( 0. .SP1 1 
Z =7 4 
NP = 1 
OC to 70 

45 X=RANDNO< C..5P1 ) a A30VE SOLAR PANEL (21 

YrWANCNCIOPO.SPZ) 

Z=Z4 
NP= 1 
GO TO 70 

EC X=P ANDN0(SP2.SP31 E UNDER SOLAR PANEL (11 

Y=RANDNO(C. .SP1 1 
Z =7 4 
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NP--1 
GO TO 7C 

55 XrRANDNOtO.»S°l ) 

Y=PANDN0(SP2.SP3) 

Z - Z4 
NP=-1 
GO TO 7 C 

er lta=randno ( c.t i.57> 

RAD=RANDNO< SP02 tSP4) 

X:R AD» COS (ETA) 

Y=RAD*STN(ETA ) 

2 - SQRT(X*X-»Y* Y)*TANTH2 +ALC2I 

70 RETURN 3 POSITION ESTABLISHES 

END 


a UNOER SOLAR PANEL 12) 


5 SECTOR 2 CONE SURFACE 
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SUBROUTINE S CP EOT { I GCC N J 

C • •*» NASA J PE 9 / 12/74 *»♦» D. EDGARS* 8I0NETICS **»• 

DIHCNSICN X t IE ) * Y ( 16 ) 

C COMMON 3L0K13 IN 0 AT A IN . EF IE-0 . ES UR r * MA I N , R NOP OS . S CPL CT . SHADE. AND CONTAM 

COMMON /BLOK 18/ SP 1 * SP 2 * S p 3 * S P 4 

C COMMON BLOK 13 CONNECTS DATA IN. EFIELC. ES'JRP.MAIM. RNDPOS AND SCPLOT. 

COMMON /BLOK 19/ SP D 1 * SPC2 * S p 1 3 

C COMMON BLOK 21 IN C HK HIT » CO NTA M. EF IELC * MA IN* RNDP OS . S H ADE * SPHIT* SCPL 

COMMON /BLOK21/ 2 1 . Z2 *Z 2 *24 , 25 

C * * • 

THIS SUBROUTINE PLOTS (CN 8.5X11 FApERJ THE MODEL AS PER THE IN°UT 
DATA TO INDICATE THE GEOMETRIC CO MPAT A8I LIT Y . THE ARGUMENT IGOON 
TELLS THE REMAINING PROGRAM TO EITHEP CONTINUE (1J OR ABORTl-ll 
BECAUSE OF UNUSABLE DATA. 

* * • 

arraV construction: 

* * • 

A- AINT (211 -1. 

X«llz se<* 

X( 2 ) = C. 

xm = o. 

X( 4 ) - SP02 
X ( 5 1 = SP4 
Xf 6>r SP03 

xm = C. 

X< 8) = C. 

X(9)r S?1 
X( 10)= SP3 
Xdll: SP 2 
X( 12) = SPC1 
X Cl 3 ) = C. 

X( 14)= c. 

X(15I=C. 

XC 16 ) =1 • 

Y C 1 ) = 22 
YC 2 ) = 22 
YC3> = 23 
YC 4 ) = 23 

Y C 5 ) = 22 
YC 6 ) = 21 

Y C7 I = 21 
YC 8 ) = 24 
YC9)= 24 
YC 10 ) = 24 

Y Cl 1 ) = 24 
YC 121= 25 

Y Cl 3 I = 25 

Y C 14) = 24 

Y CIS ) = A 
YC 161=1. 

determine: parameter usability 

IPC SP4.LE.SP02) GO TO 1000 8 PARAMETERS INCONSISTENT 

irCSP4.LE.SP03> GO TO 1CC0 S PARAMETERS INCONSISTENT 

IFCSP3.LE.SP2) GO TO 1000 8 PARAMETERS INCONSISTENT 

irCSP2.LE.SPCl) GO TO 1 CCC 8 PARAMETERS INCONSISTENT 
IFI S°2.L E.SP1 ) GO TO 1000 a PARAMETERS INCONSISTENT 
irt25.LE.Z4> GO TC 10C0 8 PARAMETERS INCONSISTENT 

ICC Z4*L E.Z 3 > GO TC 1C0C 3 PARAMETERS INCONSISTENT 

IFCZ3.LE.72) CC TC 1CC0 3 PARAMETERS INCONSISTENT 

IFCZ2.LE.21) GO TC 1CCC 3 PARAMETERS INCONSISTENT 

CC TC 16 

1C0C WRITE t 6 • 1001 ) 

1CC1 Form AT (lXf f INPUT SPACECRAFT FARAMETERS INCONSISTENT //t5Xf •RUN-- 
1AB0RTED* ) 
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I GOCN--1 
GO TO 110 
1£ CONTINUE 

* » • 

PRODUCE ELCT OF SPACECRAFT 
» » • 

2C CALL PLOTS 

30 CALL PL 0 T ( 1.C.2.G.-3) 3 ESTA3LISH ORIGIN. 

5C CALL AXIS) G..0 ..15HZ AXIS { MET ERE > » 15 » 7 . . 30 . . A,l.) a 2 
40 CALL AXISt 0. * 0. » 17HX » Y AXES ( METERS) »- 1 7 . 6 . » 0 . . 0 . . 1 . I a X.Y 
CALL LINEtX.Y. 14. 1,1.1) 8 PLOTS THE CRAFT OUTLINE. 

8C CALL SYMBOL <2 •* 5.6. .21*21 HCOME SPACECRAFT MODEL. 0.. 21) 

9C CALL SYMBOL (2.»5.2».14»25H** CINITIAL GEOM. CHECK ) **.C. .25 ) 
CALL PLOT) 12.0. 0.0.999) 

ICC IG00N=*1 
11C RETURN 
END 


AXIS DRAWING. 
AXES DRAWING. 

3 TITL 
3 TITLE 
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c 

c 

c 
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GLI3R0UT TNC S HA DC 

•••• NASA w’FL S/12/74 D.TDGARSt EXONETICS 

• • • • 

THIS routine determines whether t he particle is exposed tc the sun 

OT IS SHADED BY A PORTION CF THE SPACECRAFT. WE USE <V IA COMMON 
BLOCK l THE FIXED CONSTANTS ABOUT THE SPACECRAFT GEOMETRY OBTAINED 
DURING INFUT. 

definitions : 

(IMX.YfT) -he new particle position vector. 

C> IS E rs INPUT TC INDICATE SECTOR CF SPACECRAFT models 

(A) 1= PARTICLE If! UP’ER SECTOR. 

(B) 2= PARTICLE IN MIDDLE SECTOR. 

(C) 3- PARTICLE IN LOWER SECTOR. 

(i) ip indicates: (outfutj 

(A) lr FARTICLE IS UNDER A SOLAR PANEL. 

(B ) C= FAPTICLE IS NOT UNDER A SOLAR PANEL. 

(II ISH indicates: <outputi 

t A » lr FARTICLE IS SHADED CUTSIOE WAKE. 

(3) C- PARTICLE IS NOT SHADED. 

(C>-lr FAPTICLE IS SHADED INSIDE WAKE. 

• • * * 

REAL NEW r CS 

COMMON 3L0K02 CONNECTS CHKHIT * MAIN, RNOPOS A NO SHADE. 

COMMON /BL0KC2/ AL C 4 > 

COMMON 3L0KC3 IN CHK HIT , CONTA M. MA IN, SHADE ,S PHIT ,T RAJEC AND TRANSL. 
COMMON /3L0KC3/ NEWPCS C 3 > 

COMMON BLOK 14 CONNECTS CONTAM, EFIELD, MAIN, SHADE ANO TRAJEC. 
COMMON /8LCK14/ IF ANEL » ISH AD E 

COMMON BLCK15 IN CHK HIT . CO NTA M. EF IELD . MA IN, RNDPOS .SHADE »VNORML. 
COMMON /BLCK15/ ISECTR 

COMMON 3L0K1S CONNECTS CO NTA M, MA IN, RNDPOS ,S HADE AND VNORML. 

COMMON /BLOK 16 / NPNL 

COMMON 3L0K18 IN DAT A IN , EFIELD, ES U RF , MAI N, RNOPOS »S CPL OT , SHADE. 
COMMON /BLOK 18/ SP 1 »SP 2 ,SP 3 , SP A 

COMMON 3LOK2D IN CHKHIT, EFIELD, MAIN, RNDPOS, SHADE ANO VNORHL. 
COMMON /BL0K2FI/ T ANTH1 ,T ANTH2 , T ANT HI ,T ANT H4 

COMMON BLOK 21 IN CHK HIT , CO NTA M, EF IELD , MA IN, RNDPOS .SHADE • SPHIT, SCPL 
COMMON /BLOK 21 / 71 , Z2 » Z3 ,2 A . Z5 

COMMON BLCK23 CONNECTS EFIELD, MAIN AND SHAOE. 

COMMON /BLOK23/ XYWAKE.ZWAKE , JPCT 

EQUIVALENCE (NEWPOSIl).X) , ( NE WPOS 12) . Y! • t NEWPOS I 31 »Z ) 

EQUIVALENCE IISECTR »ISE » . t IS MDE.ISFI , ( IPANEL* IP I 
DETERMINE ISECTR 
ISE = 1 

IF t Z. GT. Z4 1 GO TO 50 
IF(Z.EQ.Z4.AND.NPNL.GE.C» GC TO 50 
IFIZ.GE.Z2) GO TO 10 
IFtZ.GE.C.I GO TO 11 
ISE=4 
GO TO 50 
1C TSE=2 
GC TO 5C 
11 ISE=3 

5C GO T0( ICC, 20 C, 300 , 3CG ) ,ISE 

** ** 

DETERMINE WHETHER PARTICLE IS ABOVE A SOLAR PANEL 


C » * * * 

ICC IFU Y. GE.SP2 .AN0.Y.LE.SP3) .AND.(X .GE.C .C.A 
1N0.X.LE.SP1) ) GO TO 11C 
I Ft ( Y.GE*C •C.AND.Y.LE.SP1) ,ANC.(X.GE.SF2.A 
1ND.X.LE.SP3) ) GC TO 11C 

C * » * » 

C MOT ABOVE/BSLCW A SOLAR PANEL. 

C * * * • 


• «* 
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115 IP = 0 

IPS ISH-C S PARTICLE IS UNSHADED. 

RETURN 

C ♦ * * * 

C IS A80VE/3EL3 W A SOLAR PANEL. 

C • ** * 

lie IR=1 

IF( ISE-21 10 6* 25C 1 26C 

C * * »* 

C PARTICLE IS IN THE MIDDLE SECTOR. 

C * * * * 

20C SPT=SP2 

C * * * * 

C DETERMINE IF PARTICLE IS SHADED IN EITHER SECTOR 2 OR 3. 

C * * ** 

2CE IFC CX.GE.SPT .AND.X .LE.SP3) . AND .CY.GE.O .C.AND . 

1Y.LE.SP1>) GO TO 110 

IFC (X.GE.C.C.AND.X.LE.SP1) .AND .CY.GE.SPT.A 
1ND • Y.LE • SP3 ) ) GO TO 110 
IPrC 

GO TO 260 

c **•* 

c particle not shaded 

c * ** * 

25C ISHrl 
255 RETURN 

C * » * 

C CHECK CIRCULAR CONSTANT 

C * * * 

26C C=X»X*Y*Y 
S=SPT*SPT 

IFtC.GT.S. AND.IP.NE.il GO TC 106 S NOT SHADED 
IF! I3E.EQ.2) GO TO 250 
ZW=SGRTIC)*TANTH4 ♦ALIM) 

IFtZ.LT.7W) GO TO 250 

ISH=-1 3 SHADED INSICE OF THE WAKE. 

GO TO 255 

C » •* * 

C PARTICLE IN SECTOR THREE. 

C * •* * 

300 SPT=SP4 

GO TO 205 
ENO 


110 


JPL Technical Memorandum 33-737 


o o oo o o o o o o r>o oo r> o o oon 


SUBROUTINE SPHIT 


SUBROUTINE SPHIT 

• *• » 

THIS ROUTINE DETERMINES WHETHER THE PARTICLE HAS HIT A SOLAR PANEL. 
DEFINITIONS: BI0NETICS«*CEE**3/13/74*» 

I'D Z IS THE NEW PARTICLE POSITION 2 COMPONENT. 

1 2) VZ IS THE NEW FARTICLE VELOCITY 2 COMPONENT. 

I 3 » Z4 IS A SPACECRAFT DIMENSION FROM INPUT PREPROCESSOR. 

I HI IS THE INDEX DENOTING: 

(A) C - DID NOT HIT PANEL. 

(cl 1 - CIC CONTACT PANEL. 

» * »* 

REAL NEWPCS .NEWVEL 

COMMON BLOK 04 CONNECTS CONTAM. MAIN. PARPOT. SPHIT AND TRAJEC. 

COMMON /BL0KC4 / COEFPP .CONST A. FRTCHC.PhIC 

COMMON BLOK 09 IN CHK HIT • CO NTA M. MA IN. S MADE »S PHIT .T RA JEC AND TRANSL. 
COMMON /BLCK09/ NEWPCS ( 3 ) 

COMMON BLOK 1C CONNECTS CONT AM . M AIN.S PHXT . T RA JEC AND TRANSL. 

COMMON /BL0K10/ NEWVEL ( 3 ) * DELTAT 

COMMON BLOK 11 In CONTAm .EFIELO. MAIN, RNDPOS .SPHIT .T RAJEC.vNORML. 
COMMON /BLOK 11 / P0SVECI3) 

COMMON 3L0K13 CONNECTS CHKHIT . CONTAM, MAIN A NC SPHIT. 

COMMON /BLCK 13/ IHITI2I 

COMMON BLOK 21 IN CHK HIT . CO NTA M. EF IELD . MA IN. RNDPOS .SH A DE » SPHI T. SCPL 
COMMON /BLOK 21/ 21 , Z2 . Z 3 .Z 4 » 25 

EQUIVALENCE ( NEWPO S ( 3 J . 2) « (NEWVEL! J) , VZ) . ( IHIT 1 1) . IH1 1 
EQUIVALENCE ( POS VE C < 3 > . ZO ) 

DZ=Z-Z4 
ADZ: ABS (DZ ) 

AVZ=ABS( VZ) 

WITHIN PANGE AND CLOSING. 

IFIDZ/VZ.LE. C.C. AND. ADZ/AVZ.LF .DELTAT) GO TO 300 

SOLAR PANEL PENETRATION. 

IFIZ.EQ.24) GO TO 300 
IF((Z0-Z4)/(Z4-Z).GT.C.) GO TO 3CC 

25C IH1=C 
RETURN 
300 IH1=1 
RETURN 
END 
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SUBROUTINE T MINf>L 

C ***» NASA J c ‘- **** 9/12/74 D. EDGARS. 8 IONET ICS **** 

C THI- ROUTINE IS A MODIFIED VERSION OF THE METE CROID IMPACT STUDY 

C SOFTWARE FREPARED IN ITS LATEST VERSION PRIOR TO THIS MODIFICATION. 

C IN DECEMBER, 1973 3Y C. BAIERLE ANC U. 5. SARENGOLTZ AFTER THE ANALYSIS 

C OF THE ''HYSICA- PRO 3L EM 3Y J. YANG. 

C ♦ » • • • • 

c 

c DECK FOR F I p ST PRESSURE FUNCTION 

C 

PARAMETER KAXFTS=5CC 

REAL I. NU 

REAL T(MAXPTS) .WCOTIMAXFTS ) .WKMAXPTS J.W2IMAXPTS) 

C COMMON 3LOK07 CONNECTS T HIN°L AND YAMG1. 

COMMON /PLC YC r / TO . V AP .FK A C C . EE » EH * ANU . E PC » R » ES » ER HC 
EOUIVALENCE ( EE.C) . IEH.HJ . (E’C. PO) . < ES.S ) . (ERHO. RHO) 

E3UIVALENCE (ANU.NUJ 

COMPUTES FIXED EXPRESSIONS FOR LATER USE 

NCELTaECC 3 MUST BE CIVISAELE BY 1C. 

ND24ANDELT-24 
ND25=NDELT-25 
NOT IN AND EL T/l 0 
AND28ANDELT-28 
AND3=NDELT-3 
AND 1C1 aN'OTEN-1 
PI a 3.14159 
ABAR : FC « S 

Da E*H»* 3/1 12.0* (1.0-NU**2) > 

BpID/(RHO*H) ) • *0 .5 
!=PI*S*A3AR 
F aS / ( 4 .0* B J *S 
Q a R/( 4. 0 * B ) * R 

COMPUTE AND STORE HDOT ARRAY 
C 

Z=ABAR/S*EXP C-R**2/S**2) 

WIDTH a TO /4. 

TEMP a Q/l WIDTH • PI) 

IF (TEMP .ST. 2.) GO TO 20 C 
N=2 

TMINrQ/(4.»PI) 

GO TO 230 

2CD IF (TEMP .GT. 48.) GO TO 210 
N=3 

GC TO 220 

21C N a: SORT ( (TEMP/12. )♦ 12.) 

22C TMIN a Q/<2 * N • PI) 

IF (THIN .GT. SORT ( P *Q I ) GO TO 230 
N = SQRT ( P *0 ) / (2.* PI) ♦ 1 
TMIN = Q / (2. • N * PI) 

23C T ( 1 ) = C.C 

WD0TC1) a: O.D 

TA = Q/l C 2*N •» 6 ) • PI) 

IF (N .GT. 3) GO TO 240 

T^AXaZ^O/PI 

CT a ( 2 . * 0/P I - TAJ/AND79 
CC TO 2 SC 

24C TM AX a (3/U?*N - 7) • PI) 

DT a ( TM AX - T A) / ANC 3 
250 T( 2) : TA - DT 
WCCTI2) : C.C 
DO 2C Ja3. ND25 
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If, °T2;p««24’A»»r 

T( J)=TA 

WA;A3AR/(RHC*H*SI* C *»2/FT2 

W3=EXPJ -P»3/o J2 > 

WCrCCS <3* T A/FT2) 

WOrSINC a*TA/FT: ) 

UOCT t J)rWA*WE» CWC + T A/r*WD) 

13 T A ; T A ♦ OT 

CONTINUE 

TF { N .ST. C> SO TO 26C 
D T = (TMAX - 2 .* 3/P I ) / 25. 

26C TA = T(N02SI*CT 

EC 25 J=NC2R. NDELT 
PT2=P**2* TA* » 2 
TC J)=TA 

WA=ABAR/(RH0*H*SI*P**2/PT2 
W3=EXPC -P*3/PT2) 

W'DrS IN ( 3* T A/PT 21 
WC=COSC 3*TA/PT2I 
WOOT( J »=WA*UB* (WC*TA/P*WD> 

TA = TA ♦ OT 
2 r CONTINUE 

C 

C COMPUTE U2 ARRAY 

C 

KEY :7 

HMlN=l.CE-8 
HMAX=1 .C 
HSTAR=HMIN 
ERMAX=l.CE-3 

DT r t TMAX - TM IN > / AND101 

TSAVE=T(2) a FOR VAP CALCULATION BELOW. 

TLO = Tf 2 ) 

THI ; THIN 

DO 30 J = l, NDTEN 

CALL ROMBS (TLO.THI »TT »FOFTT *HST AR » HMIN *HMAX .ERMAX • ANStK * KEY | 

60 FC=0. 0 

FC = -2. / (TO • (1 ♦ ( T HI-T T) /T C) • *2 J 
70 TW=TT 

CALL SLUP (TW.FCC.FP.T.WOOT, KDELT.2I 

FOFTT=FC*FCC 

CALL R0M2 

IF CK.E3.1) 60 TO 60 

W2» J> =ANS*2 .59E-6 S TO CONVERT TO KILOGEE FROM ENGLIAH. 

THIrTWI* OT 
SC CONTINUE 
C 

C COMPUTE W2MAX 

C 

ice u:max=c.c 

EC 11C J=1.NET"N 

IF CW2(J> .LE. W2H.AX) CC TC 1 1 C 
U2MAX = W(J» 

JACK=J 

110 CONTINUE 
C 

C THIS CALCULATES 1 VALUE OF VA° AT TIME WHEN THE ACCELERATION IS 

C AT THE MAXlMUN VALUE. 

C 

TLO=TSAVE 
A J = JACK 

THIrTMlN-M AJ-1 ,)*DT 

CALL ROMSS (TLO. THI. TT , FO FTT. HS TAR , H MIM, HMAX , ERMAX. ANS. K. KEY J 
FS : 1. / (Cl ♦ (THI-TT)/TC)»*2) 
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SUBROUTINE THINPL 


4T TW=TT 

CALL SLUP (TV»FSS»FP.T»WD0T»NDELT»2> 

FOFTT=FS*FSS 
CALL R0M2 

VAPrANS* .0254 51 CONVERT TO M/SEC 

CHECK SIGN OF VAP 

IF< VAP.GE.C.) GO TO 525 
GO TO 12C 

COMPUTE VI ARRAY WHEN TRIAL IS POSITIVE. 

£25 DT=(TMAX-TMINJ/AN01C1 
TLO = TSAVE 
THI = thin 
DC 55C J = 1»NDTEN 

CALL ROMBS (TLC »T HI *TT » FOFTT .HST AR • VMIN»HMAX .£ RMAX • ANS»K » KEY ) 
530 FS=0.0 

FS = 1. / (Cl ♦ CTHI-TT>/TC>»*2) 

540 TW=TT 

CALL SLUP CTV»FSS»FP*T .WDOT • NDELTt 2 ) 

FOFTT=FS*FSS 
CALL R0M2 

IF (K.EQ.l) GO TO 530 

WHJ)=ANS*.C254 9 CONVERT VELOCITY TO M/SEC 

THI=THI*DT 
55C CONTINUE 
VAP=0. 

W2HAX=0. 

DO 60C J = 1*NDTEN 
IFIU1CJ».GE.C. J GO TO 500 
IF( W2( J ) .LE.W2MAX) GO TO 600 
W2MAX=U2t J) 

VAP=W1< J) 

60C CONTINUE 
12C PK ACC = W2MAX 
RETURN 
END 
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c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

w 

c 

c 

c 

c 

c 

c 

c 


c 

c 

c 

c 

c 


c 

c 

c 


c 

c 

c 


c 

c 

c 


c 

c 

c 


SUBROUTINE TRAJEC 

♦ NASA v’PL •••• 9 /ID/ 74 ***• D. EDGARS* BIONCTICS **** 

***•• ***** ***** 

SUBROUTINE TRAJEC CALCULATES THE INCREMENT IN THE PARTICLE TRAJ- 
ECTORY* THIS REQUIRES INFORMATION ABOUT THE CURRENT PARTICLE 
CHARGE AT THE LOCATION TO BE MOVED FROM* THE SUN/SHADE CONDITION* 

THE CURRENT VELOCITY VECTOR /ND THE ELECTRIC FIELD FROM THE CHARGED 
SPACECRAFT. CAL*. S TO SUBROUTINES EFIELD AND PARPC T PROVIDE THE 

field and charge information, this subroutine calculates a suitable 

TIME INTERVAL c OR THE STANDARD KINEMATICAL EQUATIONS OF MOTION 3Y 
PRECALCULATING A DELTA T FOR A VELOCITY DEPENDENT* FORCE -FIELD 
DEPENDENT AND A DRASTIC ELECT 0 IC POTENTIAL CHANGE DOMINANT CASE AND 
USES THE SMALLEST OF THE THREE. BY THIS METHOD WE MAINTAIN A FAIRLY 
UNIFORM FOUR CENTIMETER INCREMENT SIZE EXCEPT WHEN THE ELECTRIC 
POTENTIAL DIFFERENCE WOULD BE TOC LARGE OVER THE DISTANCE SPECIFIED 
(4 CM) AND THUS* IN ESSENCE* WE SMOOTH THE SUN/SHADE SOUNDRY 

traversal. 


***** ••••• ••••* 

DIMENSION RPCSC31.RVELI3I » OELTAR C 3) . C POS C 3 ) . DV EL C 3 > .EC3I 
COMMON BL0KC4 CONNECTS CONTAM. MAIN* FAPPCT* SPHIT AND TRAJEC. 
COMMON /BLCKC4/ CO EFRP .CONSTA » PRT CH3 • PHIC 

COMMON BLCKOP IN CHXHIT * CO NT iK*MAIN*SHACE* SPHIT* TRAJEC AND TRANSL. 
COMMON / 9LCK 0 3/ NEWPCSC3) 

COMMON BLOK1C CONNECTS CONT A M* MA IN *S PHIT »T RA JEC AND TRANSL. 

COMMON / BL OK 1 0/ NEWVEL13). DELTAT 

COMMON BL0K11 IN CONT AM. EFIE LO *M AIN. RNCPOS .SPHIT » TRAJEC* VNORML . 
COMMON / BL OK 1 1/ POSVEC 13) 

COMMON BLCK12 CONNECTS CONTAM. MAIN AND TRAJEC. 

COMMON /BLOK 1 2/ SMASS .VELVECC3) .ILCST 

COMMON BL0K14 CONNECTS CONTAM. EFIELD. MAIN. SHADE AND TRAJEC. 
COMMON / BLOK 1 4/ IP AN EL ♦ ISHADE 

COMMON BLCK2E CONNEC T S CONTAM. EFIELC. MAIN AND TRAJEC. 

COMMON /BLOK 2 5/ E VEC C 3) 

COMMON BL0K27 CONNECTS PARPOT AND TRAJEC 
COMMON /3L0K27/ PHICN6 

EQUIVALENCE CISHADE.ISH) . « EV EC 111 . E Cl >1 
EQUIVALENCE CRPOSIl) .POSVEC 11 )> . (R VEH 1) ♦ V ELV EC 111 ) 

EQUIVALENCE CDPOS Cl I .NEWPOS CD). COVELC 1 I .NEWVELC 1 1 J 
CALL EFIELD 

CALCULATION OP THE POTENTIAL DEPENDENT TIME INTERVAL. 

CALL PARPOT 
OELT3=10. 

IFC ABSC PHICNO.LE.O.l) 50 TO 2 
DELT3=1./ABS CPHICNG) 

2 TFRP=0. a IN SHAOE 

IFCISH.EO.O) TPRP = COEFRP a IN SUN 

VELOCITY DEPENDENT DELTA T. 

CELT 1 =C.C4/ CSQRTCRVELC1 J*RVELC1MRVELC2» *RVEL C2 MRVELC3J*RVELC3)> 

1*1. E-6) 

FORCE DEPENDENT DELTA T 

0ELT2=SQRTC3.08*GN ASS/ CA3S CTF9PI «-ABS C PRT CHG )*SQRTIEC1)*EI1I* 
1EC2»*EI2)+E »3»*EC31 M1.E-2C) ) 

THE MINIMUM DELTA T OF THE THREE CHOICES 
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SUBROUTINE TRAJEC 


DELT -AM INI CDELT1 t DELT2 ) 

DELTATrAMlNlf CELT t DELT3) 

5 CTSOTM=DELTAT* DELTAT/C2.*GMASS ) 

K3 -0 

1C DO 15 Krl f 3 

IF ( K. E Q • 3) K3rl 

15 CELTAP(KJ=DTS0TM*(TFPP*K3+ PFTCHG*E CKI 1+RVELCK JOEL TAT 
DO 4C K-l # 3 

IFCCELTAR(K) .GT.C.04) GO TO 5C 
4C CONTINUE 
GO TO 31 

50 0CLTAT=DELTAT/2. 

GO TO 5 

31 00 35 K-l t 3 

DFCS CK JrRPOS CK ) ♦ DELTAR { K) 

35 DVELf K)=2.*0ELTAR CK I / DELT AT-* VEH K) 

DELTAGrCONSTA*PHICNG*DELTAT 9 CHANGE IN PARTICLE CHARGE. 

PHIQ=PHIO*PHICNG*DELTAT 9 NEW PARTICLE POTENTIAL. 

PRTCHG=PRTCHG*CELTAQ 

RETURN 

END 
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SUBROUTINE TP AN SL 

C NASA JPL •••• 9/12/74 0#ECCARSt BIONETICS 

C THIS ROUTINE TRANSLATES THE p ARTICLE POSITION A NC VELOCITY BACK INTO 

C THE FIRST OCTANT BECAUSE OF SYKKET RY TO EASE TWE CALCULATIONS 

C IXtYtZ) IS THE POSITION WHICH REQUIRES TRANSLATION - (NEW COORDINATES) 

C CVXiVYIrTHE VELOCITY WHICH REQUIRES TRANSLATION - CNEW VELOCITY) 

C * • • • • * 

REAL NEVPCSt NEWVEL 

C COMMON 3LCKC9 IN CHK HIT • CCNTA M. MA IN. S HADE t S PHIT tT PA JEC AND TRANSL# 

COMMON /BLOKS' 5 / NEUFCS C 3 ) 

C COMMON 3L0K1C CONNECTS CONTAKt MAINtSPHlT .TRA JEC AND TRANSL# 

COMMON /eLOKIC/ NEWVEL C 2 ) • OELTAT 

EQUIVALENCE l XtNEWPOS Cl) ) t(Y*NEWP0S(2>) t i NEWVEL Cl >tVX It (NEWVEL<2)t 
IVY) 

IFC X#LT.C.C.AN0.Y#GT#C.C) 00 ^0 ICC 
IFCX.LT.r.C.AND. Y#LT#C.D GC TC 2DC 
IF ( X# GT#C#n#AND#Y#LT#C#D 30 TO 3CC 
RETURN £ DCES NOT ALTER PARAMETERS 
ICC XP=Y 
YPr-X 
VXP-VY 
V YP - - V X 
GO TO 4CC 
2CC Xr-X 
Y--Y 
VXr-VX 
VY=-VY 

return 

3CC XP --Y 
YP-X 
VX°-- VY 
VYP= VX 
HOC X=XP 
Y=YP 
VX=VXP 
VYzVYP 
RETURN 
END 
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c 

c 

c 

c 

c 

c 

c 


SUBROUTINE VNORML 

NASA JPL ♦*** 9/12/74 D. EDGARS* EIONETICS ***• 

* * * * * * 

THIS ROUTINE CALCULATES THE NORMAL UNIT VECTOR AT THE SURFACE IMPACT ZONE 
THE VECTOR TS ALONG (Z-AXIS) FOP SOLAR PANEL CONSIDERATION. 

(IT THEN ADDS A ID PERCENT BIAS TO ONE OF THE VELOCITY COMPONENTS.) 

NP-*1 ... A3CVE PANEL: -1 ... 9EL0W PANEL. 

* * * 

DIMENSION VELCn * TANT<4) 

COMMON BLOK11 IN CONT AM * EFIE ID t M AIN t RND PCS t S PHIT * TR A JE C . V NO RML . 

COMMON / 3LOK1 1/ P0SVEC<3) 

COMMON BLCK15 IN CHKHIT .CO NT /M.EFIELC* MAIN. RNDPOS .SHADE. VNORML. 

COMMON / BLOK1 5/ ISECTR 

COMMON BLOK1S CONNECTS CC NT AM t M AIN . RNCPOS .SHADE AND VNORML. 

COMMON /SLOK1S/ NPNL 

COMMON BL0K16 IN OAT AIN *EF IE LD *ESU RF • MAIN • RNCPCS . SCPLO T* SHADE . 

COMMON / 3LOK 1 9/ SPltSP 2 * SP3 *SP4 

COMMON BLCK2C IN CHKHIT* EFIELOt MAlNt RNOPOS* SHADE AND VNORML. 

COMMON /3L0X20/ T ANTHI . TANTH2 . T ANTH3 tT ANT H4 

COMMON BL0K2I IN CHKHIT.COnT AMtEFIELC. MAIN* RNDPOS .SHADE. SPHI Tt ’sCPL 

COMMON / BLOK 2 1/ Z1 . Z2 . 23 t Z4 » Z 5 

COMMON BLCK22 CONNECTS MAIN ANO VNORML. 

COMMON / 8L0K2 2/ VELNRM I 3 ) 

EQUIVALENCE CNPNL *NP > . C VELNRM < 1 > tV EL C 1 >) 

EQUIVALENCE i PC SVEC 1 1 ) * X) .(PDSVECI2) * Yl 
EQUIVALENCE CPOSVECf 3J »Z> 

TANK 1 ) rTANTHl 
TANT(2)=TANTH2 
TANTC 31 r TANTH3 
T ANT ( 4 I rT ANTH4 
XYOI$T~SQRTlX*X*Y*Y) 

IFCNPNL.NE.O GO TC IOC 

IS -4-ISECTR 

FXr-X* T ANT (IS1/XYDIST. 

FY ~-Y*T ANT ( ISl/XYDIST 
FZr+1. 

DE NOMrSQRT < FX *FX*FY*FY*FZ*FZI 
VEL<ll=FX/DENOM 
VEL(2)=FY/DEN0M 
VELC3>=FZ/DENOM 
GO TO 300 
ICC VELfllrC. 

VEL<2UQ. 

RNzNP 

VELC3URN a ♦ OR - Z DIRECTION AS PER NP 

3CC CONTINUE 

NO=RANQNO<l.t 3.99) a ADO TEN PERCENT TO ONE OF THE COMPONENTS. 

GO T0(31C*32O. 3301. NO 
31C VELf D-VELI11+.1 
GO TO 400 

320 VEL«21rVEL<2) >.l 
GO TO 400 

33C VELI3UVELI3I *.9 
4CC RETURN 
END 
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SUBROUTINE YANG 1 


SUBROUTINE YANG1 a THE SURFACE IMPACT ANALYSIS 

C **** NASA JPL »*«• S/12/74 **•* 0. EDGARS* BIONETICS 

C THIS ROUTINE IS A MODIFIED VERSION OF THE METEOROID IMPACT STUOY 

C SOFTWARE PREPARED IN ITS LATEST VE RS I CN , PRIO R TO THIS MODIFICATION. 

C IN DECEMBER. 1373 BY C. BAUER LE AND J. 3. 3 ARENGOLTZ AFTER THE ANALYSIS 

C CF THE FHYSICAL PROBLEM BY J. YANG. 

c • * • • • • 

c 

c 

C DR*- - BASIC DATA ARRAY FOR RADII 

C RR — RADIUS ARRAY FOR INITIAL POINTS 

C FR— FINAL RADIUS ARRAY 

C RPKACC— INITIAL ARRAY FOR PEAK ACCELERATION AT VALUES IN RR 

C FPKACC — FINAL ARRAY FOR PEAK ACCELERATIONS FOR VALUES IN FR 

C RVAP— ARRAY FOR VELOCITIES AT PEAK ACCELERATIONS AT RR RAOII 

C FVAP — FINAL ARRAY FOR VEL. AT PEAK ACCELERATIONS AT FR RADII 

f* 

w 

C ••••••*••••• 

c 

c 

real nu 

DIMENSION RR (14 I * RPKACC(14I» RVAP(14I» 0RC14 ) 

DIMENSION FR (35) » RMID<2>. AHID«2>. 

2 VMID ( 2 ) 

C COMMON BLOK 05 CONNECTS OATAXN* MAIN. RELEAS A NO YANG1 

COMMON /BL0KC5/ KO .PM ASS .PRH C.RS (35 J . FV AP ( 35 ). FPKACC (35 )» RHO 
C COMMON BLOK 06 CONNECTS OATAXN* MAIN AND YANG1. 

COMMON /BL0K06/ E»H*NU»PO»ISKIP 
C COMMON BL0K07 CONNECTS THINPL ANO YANG1. 

COMMON /8L0K07/ TO*VAP*FKACC»EE*EH*ANU*EPO*R*ES»ERHO 
C COMMON BLOK 17 CONNECTS CATAIN* MAIN ANO YANG1. 

COMMON /BL0K17/ PVEL»ACCNIN 

equzvalence(fr(1)«rs( 1)) a yangi 

EQUIVALENCE (KO.IMAXI) 

DATA DR / .01* .02* .03* .04* .06* .08* .1* .2* .4* .7# l.t 
2 0.*C.*0. / 8 THE OLD COOING HAD INTEGER ZEROES. 

10D0 FORMAT! 1H0. IX. ’THE PLATE HAS PUNCTURED*) 

ANU=NU 

N0ELT=50O 

ISKIP=C 

C 

WRITE! 6*4701 PMASS .PVEL.PRHO .E.H.NU.RHO.PO 
470 FORMAT! 1H0* *T HE FOLLOWING VALUES ARE RESPECTIVILV**/*1X**-PNASS — P 
IV EL — PR HO — E — H — NU--RHO — PO *./ . IX . 8 (1PG9.4 *2X )) 

DO 7 1=1. 14. 1 

RR(I) = DR(I) 

7 CONTINUE 

PI = 3.14153 

CALCULATE PR-P ARTICLE RADIUS* TO* AND S 
PR = CBRT ( (3 . • PMASS) / (4. • PI • PRHOI I 

DETERMINE WHETEHER CONDITIONS EXIST FOR PLATE PUNCTURE BV METEOROID. 

P=CBRT(12.*PMASS*PVEL»PVEL/f PI*P0))«.5 

IF(P.LE.H) GO TO B 
C PLATE PUNCTURED. 

TC = 8./3.*PR»(P0/(2.»PRH0))*i(l./6.)*PVEL**(-4./3.)»ll.«- 
1 3 •* R HO * H/P R* CBRT ( ( P VEL/2 . ) » • 4/ ( PRHO»FC*PC 1 ) ) 

WR ITE ( 6 * 1C00 ) 

GC TO S 

PLATE NOT PUNCTURED. 
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SUBROUTINE YANG 1 


e tc = e. * pr • (x. ♦ o. » phc * fvel • fvel ) / pci 

2 * 1C.» • ( Ai.0GI0(?r/(2 .•PRJIOI ) /E. - I 4./ 3.1* AL0G19(PV£LI 1/3. 

5 CONTINUE 

S = (PR/3.) * C9R T ( ( 2 • * RHC • PVFL • PVEL) / PC) 

J = 5 


ACCMlN IE THE MlNlV'JM ACCELT SATICN CF INTEREST 

VARIABLES CHANGED F^CM METIC TO ENCLISP SYSTEM 

EE=E»1 .45E-4 
EH=H»39. 37 
ERHO:RhO*9 .3S6E-3 
EP0=P0*1.45E-4 
ES-S*23.3 7 

VALUES FOR PEAK ACCELERATION AND VELOCITY AT PEAK ACCELERATION ARE 
CALCULATED FOR INITIAL VALUES OF RADII. AT LEAST 2 RADII MUST 
BE USEC. AFTER THAT IF THE LARGEST VALUE FOR A PEAK ACCELERATION 
C HAS BEEN PASSED AND THE PEAK ACCELERATION ,JS LESS THAN ACC MIN* 

C CALCUALTIONS STOP. 

C 

DC 1C J = l* 11 

R = RR(J) 

R = R • 39.37 
CALL THINPL 
P = R / 39.37 
RPKACCtJ) = PKACC 
R VAP ( J ) = VAP 

IF ( J .LT. 2) GO TO 10 
IF t RPK ACC ( J-l I ,L T . PKACC I GO TC 1C 

IF (PKACC .GT. AC CM IN I GO TO 10 

JJJ : J 
GO TO 15 
1C CONTINUE 

JJJ = J 
15 CONTINUE 
11 CONTINUE 

ACCMAX = RPKACCtJJJl 
KKK = JJJ 
C 

C SEARCH FOR THE MAXIMUM VALUE OF PEAK ACCELERATION IF AVAILABLE. 

C IF IT IS NOT » THREE MORE VALUES OF RAOII ARE CHOSEN TO FIND IT. 

C ONLY THE ONES NEEOED ARE USED. WHEN MAXIMUM IS FOUNO NO MORE 

* C ARE CALCULATED. 

DO 2C J= 2. JJJ 
I = JJJ ♦ 1 - J 

IF (RPK ACC (I I .LT. ACCMAX I GO TO 35 

ACCMAX = RPKACC(I) 

KKK r I 
CC CONTINUE 

KKK : 2 

OT ST : ( R 7 (1J - S*l.:>/3. 

CO 3C I rl, 3 

R = RR( 1 I - D 1ST 

R = R • 3F.37 

CALL THINPL 

R = R / 39.37 

JJJ = JJJ ♦ 1 

00 25 J = JJJ. 2. -1 

RR ( J I = PR(J-l) 

RPK ACC ( Jl : PP K A CC ( J-H 
RVAP(J) = RVAP(J-l) 
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25 


C 

c 

c 

c 


32 


35 

36 

C 

c 

c 

c 

c 

c 

c 

c 

c 

c 


4 C 


45 


CONTINUE 
RRtll = r R 
RPKACCC1) =PKACC 

RVAPfll 3 VA© |( 

IF f RPK ACC (1 I .LT. ACCMAXI CG TC 35 
ACCMAX 3 PKACC 
3C CONTINUE 

IF THE MAXIMUM PEAK ACCELERATION IS NCT FOUND THE FINAL ARRAYS 
FORMED APE PRINTED AND CALCULATIONS CEASE FOR THIS CASE. 

IMAXI=JJJ 

DC 32 Irl, JJJ. 1 
FRI I|3R<M I) 

FPKACC(I) rRPKACC(I) 

FVAPC I):RVAP( II 
CONTINUE 
GO TO T9 
CONTINUE 
CONTINUE 
LLL 3 C 

A FIBONACCI SEARCH IS DCNE TC FIND THE MAXIMUM PEAK ACCELERATION. 
THE TWO VALUES BRACKETING THE MAXIMUM PREVIOUSLY FOUND ARE USED AS 

the starting end points, ruffer — largest value OF R» AuPPER and 

V’JPPER CONTAIN THE CORRESPOND INS PEAK ACCELERATION AND VELOCITY. 
RLOWEP--LCUE* VALUE OF R. R MID ( 1 >- -LESSER CENETER VALUE OF R. 
RMIOC 2) — THE GREATER CENTER VALUE OF R. 

WHEN A VALUE OF RUPPER IS CISCARDEO FROM THE SEARCH IT IS KEPT 
IN THE FR ARRAY AS PART OF THE FINAL RESULTS. 

RUPPER 3 RRCKKK *1 1 
AUPPER3RPKACCCKKK*1> 

VUPPER = RVAP (KKK*1I 

RLOWER = RRCKKK-ll 

RDIFF = RUPPER - RLOWER 

R r RLOWER ♦ .381967 • RDIFF 

R r R • 39.37 

CALL THINFL 

R = R / 39.37 

RMIDI1I 3 R 

AMID! II 3 PKACC 

VMIOC1I 3 VAP 

III 3 1 
LL 3 2 

R 3 RUPPEP - .331957 • RDIFF 

CONTINUE 

III 3 III ♦ 1 

IF (III .ST. 11) GO TO 6 C 

R 3 R • 39.37 

CALL THINFL 

R 3 R / 33.37 

RMID(LL) r R 

AMIOtLL) 3 PKACC 

VM ID (LL 1 3 VAP 

IF ( AM ID ( 1 ) .GT. A* ID < 2 ) I 3C TO 50 
RLOWER 3 RMIDC1I 
RDIFF1 3 RUPPER - RMIDC2I 
RDIFF2 3 RMI0C2I - ©LOWER 
IF ( RDIFF1 .GT. RDIFF2I 30 TO 45 
R 3 RLOWER ♦ RCIFF1 
LL 3 l 
GO TO 4C 

RMIDC1I = RMIDC2I 
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SUBROUTINE YANG 1 


AMI0I1I - AMID (2 1 
VM IOt 1 ) = VMID<2) 

R = RUPPER - RDIFF2 
LL = 2 
GO TO 40 

50 LLL = LLL ♦ 1 

IF (LLL .EG. 1J GO TO 54 
DO 52 1= LLL. 2. -1 

FR ( I) = FR(I-1> 

FPKACC(I) = FPK ACC ( 1-1 ) 

FVAP(I) : FVAP (I-lf 

52 CONTINUE 

54 FR ( 1 ) = RUPPER 
FPKACCtl) = AUP°ER 
FVAP (1 ) = VUFPER 
RUPPER = RMIDC2) 

AUPPER = AMIC { 2 ) 

VUPPER = VMI0I2J 

RDIFF1 = RMIcm - RLOWER 
RDIFF2 = RUPPER - RMIDU) 

IF (RDIFF1 .GT. RDIFF2 I GO TO 55 
R r RUPPER - RDIFF1 
LL = 2 
GO TO 40 

55 RMID (2 ) = RMIO(l) 

AM I0( 2 ) = AMID(1» 

VMID(2) = VMID(l) 

R = RLOWER ♦ R0IFF2 
LL = 1 
GO TO 40 

"1 IF (LLL .EG. 01 GO TO G7 
DO 65 1= 1. LLL. 1 

J = I ♦ 10 
FRIJJ = FR( II 
FPK ACC C J) = FPKACCCI) 

FV API J ) = FVAPII) 

4.5 CONTINUE 

RD IFF : tRUPPER - S)*.l 

THE FINAL ARRAY INCLUDES TEN E3UALLY SPACED VALUES FOR RADIUS 

between the maximum feak acceleration and the crater wall. 

DC 7C 1= 1. 1C. 1 
FR ( I) : F?( 11 ) - (11 - I) • RDIFF 
FPK ACC ( II = FPKACC(ll) • FR (11 1 / F R (I I 
FVAP(I) = FVAPtll) * FR(ll) / FO(I) 

CONTINUE 

VALUES FROM THE INITIAL ARRAYS ARE ACCEC TC TPE FINAL ARRAYS. 

LL = KKK ♦ 1 
CO 75 1= LL. JJJ. 1 

J - 9 ♦ LLL ♦ I - KKK 

FR ( J I : RRII) 

FFKACCIJ) = PPKACCdl 
FVAP(J) = RVAP(I) 

75 CONTINUE 

M = 9 ♦ LLL ♦ JJJ - KKK 

IMAXIrM 
79 CONTINUE 
81 WRITE l 6. 3C0 > 

300 FORMAT(IHC) 

DO 250 Mlrl.KC 
IFIFPKACC(Ml) .LE.C.1 ISKIP=1 


122 


JPL. Technical Memorandum 33-737 


SUBROUTINE YANG1 


25C WRlTE<e»Z75) M 1, FR !M1» *M1» FPKACC1M1 I .Hi t FVAP IMl I 
275 FORMAT! 1X»*FR !• » I3t * > = • « G9.*» 2X * ‘FPKACCI * »I3» • » = *# 69. R* 2X» 
1’FVAPf *»I3»*lr»,G9.AI 
I0UMMY=0 
00 260 IY=1»K0 

WRITE! 7t 16011 IY» IOUMMY.FR ! IY » t FPKACCI IY) .FVAPHYI 
260 CONTINUE 

1601 FORMAT! 2I3.1P3E15. 5) 

RETURN 9 WAS 60 TO 5 
END 
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The general model processors have three basic forms: (1) the fully 
deterministic complete driver with intermediate data output (FULDET), 

(2) the driver (FULBPS) that bypasses the meteoroid model/ surface response 
model and operates on the intermediate data output generated in a type (1) 
execution, and (3) the driver for the pyrotechnic simulation (PYROM). Pre- 
sented below are the computer control card run streams with the proper 
sequence for the execution of a case study. In each example, the data cards 
are preceded by a line indicating the format. The basic assumption here is 
that the model exists on a magnetic tape compatible with a Univac 1108 major 
computer system. On the tape are three files corresponding to FULDET, 
FULBPS, and PYROM processors respectively. 

Example case studies are given below: 

(1) The fully deterministic model, FULDET, allows the sequential 
selection of meteoroid mass and velocity groups with the option 
of setting the range of dust-grain diameters and the number of 
random impact positions. The data card used for this run is of 
the form 

MML MMU MVL MVU IGML IGMU N2 IHICUP 

where MML and MMU are the meteoroid -mass index limits (may 
use integers from 1 to 10, as in Table 3, MVL and MVU are the 
meteoroid -velocity index limits (may use integers from 1 to 3, 
as in Table 3), IGML and IGMU are the ejecta -diameter limits 
(may use integers from 1 to 10 corresponding to 10 to 100 pm), 

N2 is the number of positions to be randomly selected per each 
meteoroid mass -velocity group, and IHICUP increments the 1108 
system random number generator (because the computer produces 
the same set of random numbers for each case, unless manually 
incremented on to new numbers). A set of acceleration- and 
velocity vs. range data cards is output for each mass -velocity 
group (for use with FULBPS). Table D-l provides a sample run- 
stream where the contents of file one of the magnetic tape are 
copied, the element DA TAIN is updated for a particular mission 
when appropriate (cards from the first @ DELETE to the @ PREP 
card may be omitted when DAT AIN is acceptable), and execution 
commences. This particular run is for meteoroid-mass groups 
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2 through 10, meteoroid-velocity group 3, grain-diameter of 
10 \xm only, and 1 random position* This would be typical for a 
run that is meant to produce only the acceleration and velocity 
vs* range data cards. 

(2) Similarly, Table D-2 depicts a FULBFS runstream for a case 
study involving meteoroid -mass group 3, velocity group 1, and 
all 10 ejecta diameter possibilities for 100 random positions. 

The three types of data cards (format indicated in Table D*-2) arc! 

N2 KPLOT 

with N2 (integer) positions, KPLOT index indicating whether to 
plot (via Calcomp) the spacecraft outline for the integer 1 or to 
skip the plot when 0 is entered, 

IHICUP KO IGML IGMU PVEL PMASS IMV IML SMI 

where KO is the dimension of the following data array, IGML and 
IGMU are the grain diameter limits (1 to 10), PVEL is the mete- 
oroid velocity, PMASS is the meteoroid mass, IMV and IMM are 
the meteoroid velocity and mass indices (as per Table 3) and 
SMI is the corresponding value related to the flux and velocity of 
this particular M-V group by 

SMI = FLUX (IMM, IMV)*12/PVEL, (Dl) 

and finally the acceleration and velocity vs. range array where 
the data is in the form of 

I R(I) A(I) V (I) 

with KO entries. 

(3) Lastly, the pyrotechnic event runstream appears in Table D-3. 
The only essential difference between this runstream and the 
one for FULDET is the data card preceding the acceleration and 
velocity vs. range array. The form of this card is 

IHICUP KO IGML IGMU IPYRO NPYRO PRMIN PRMAX ANAME 

where IPYRO indicates the identifying number assigned to this 
type of device, NPYRO indicates how many pyrotechnic devices 
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of this type are on the spacecraft, PRMIN and PRMAX are the 
localizing parameters defined in Appendix B and ANAME is an 
18 -space alphanumeric spacecraft zone name used for printout 
purposes (see Table 2b for example). 
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TABLE D1. SAMPLE FULDET EXECUTION RUNSTREAM 


3RUN FULDET » ETC. 

SASCtT TAPE. *T»XXXXX 
3REWTND TAPE. 

8CCP Y »£ TAPE.tTPFS. 

3FREE TAPE. 

3DELETE t A .FULAES 
aOELETE* SR .DATAlN 
aFOR.lS .DATAlN, .DATAIN 
(THE NEW DATAIN DECK) 

3PREP 

3MAPtTS .FULMAPf.FULABS 
IN .FULDET 
IN .DATAIN 
lie LIB*JPLS. 

LIB LI 5 *PLCT $. 
aXGT .FULABS 
• * FORMATt 3(13* 3X) ) • * 

2 1C 331 111 

3FIN 
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TABLE 02. SAMPLE FULBPS EXECUTION RUNSTREAM 


•RUN FULBPS »ETC. 

•A SO • T TAPE. .T.XXXXX 
•REMIND TAPE. 

•NOVE TAPE. *1 
•COPY *6 TAPE.*TPFS. 

•DELETE* SR .DAT AIN 
•DELETEvA .FULBPA 
•F OR*IS .DATAIN*. DATAIN 

(THE NEW DATAIN DECK! 
•PREP 

•MAP *IS .FULBPM*. FULBPA 

IN .FULBPS 


IN 

•DATAIN 



LIB 

LIB* JPL$. 



LIB 

LIB*PLOT$. 



• XQT 

•FULBPA 



• • 

FORMAT ( 2110 1 • • 




100 n 



• • 

FORMAT! 4(13*2X1 *10X 

*2E12.5*1X*2CI2 

*2X!*E9 .41 • 

1 

19 1 10 

• .14453+05 .10000- 

• • 

F0RMAT(I3*3X*1P3E15 

.51 * * 


1 

7. 56170- C 5 

1. 674 86 +G1 

-1 .43073-01 

2 

8. 73974-04 

1 • 44 510 + 00 

-1.23788-02 

3 

1.67233-03 

7.57313-01 

-6. 46925-03 

4 

2.47C69-C3 

5. 12601-01 

-4.37883-03 

5 

3.26904-03 

3.87415-01 

-3. 3C 94 5-03 

e 

4.00740-03 

3.11373-01 

-2.65986-03 

7 

4. 86576-03 

2.60284-01 

-2.22344-03 

8 

5.66411-03 

2. 2359 7- Cl 

-1.91CC5-03 

9 

6.46247-03 

1.35974-01 

-1 .6743 3-03 

10 

7.26083-03 

1 . 74426-01 

-1 .49CC1-03 

11 

3.05919-03 

1.57147-01 

-1 .3424 1-G3 

12 

8. 20090- 03 

1.54535-01 

-1 .32672-03 

13 

3.42989-0 3 

1.5}27 37-01 

-1 .31693-03 

14 

9*02959-03 

1 . 4 0 6C 3- 0 1 

-1. 21274-03 

15 

1. rGCCC-C2 

1 .24262-01 

-1 .06287-03 

16 

2. 00000-02 

6.59525-02 

-6.16022-04 

17 

3.rCCCG-C2 

4.44363-02 

-3 .83329-04 

18 

4.00000-02 

3.44277-02 

-3.00458-04 

19 

6. CCCCC-C2 

2.33315-02 

-2.04667-04 
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TABLE D3. SAMPLE PY ROM EXECUTION RUNSTREAM 


aRUN PYROM»ETC. 

3AS6*T TAPE. »T .XXXXX 

3REWIND TAPE. 

3MOVE TAPE. # 2 
aCOPY*G tape.*tpf$. 
aFREE TAPE. 

BDELETE * A .PYROMA 
iOELETE»SR .DAT AIN 
BPORflS .DATA IN*. DATA IN 

(THE NEW DATAIN DECK) 
apREp 

aMAP*IS .PYROMS*. PYROMA 


IN 

.PYROM 


IN 

.DATAIN 


LIB 

LIB* JPL* • 


LIB 

LIB*PLOT S. 


* • 

FORMAT ( 2I1C 1 • * 



100 0 


• * 

FORMAT ( 6d3»2X)t2El2.5» 

3 AG) • * 

002 

06 Q01 10 001 004 

.22*00 

• • 

FORM AT Cl3f3X»lP3E 15.51 

* * 

1 

5.2-04 

1.00*02 

2 

8.0-C2 

2.80*C1 

3 

2.5-C1 

3.0C*0C 

4 

5.0-01 

S.0C-01 

5 

T.5-C1 

4.CC-C1 

6 

1.00*00 • 

1.00-01 


8FIN 


• 54 ♦DOS £CT OR ONE CONE 


-1. 55+01 
-4 • EC*CC 
-2.7C-C1 
-6.3C-C2 
-3.SC-C2 
-2.70-CZ 
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